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FOREWORD 

This  work  reflects  ongoing  research  that  has  been  sponsored  within  the  Combat  Systems 
Department,  Naval  Surface  Warfare  Center,  Dahlgren  Division,  that  responds  to  the  need 
for  improved  techniques  in  dealing  with  the  situation  where  many  threats  are  attacking  a 
battle  group.  In  such  a  situation,  each  threat  has  a  large  number  of  potential  targets  from 
which  to  choose.  Prioritization  of  the  threats  within  a  raid  that  involves  a  large  number 
of  Anti-Ship  Cruise  Missiles  (ASCMs)  presents  a  problem  to  a  global  scheduler  in  terms  of 
both  selection  of  shooting  platforms  and  determining  the  queue  for  engagement  scheduling. 
This  problem  becomes  particularly  challenging  when  placed  in  the  context  of  the  network 
centric  paradigm  envisioned  for  the  future  of  the  surface  navy,  and  represents  a  significant 
potential  problem  that  needs  to  be  addressed. 

This  work  not  only  addresses  the  problem  by  discussing  its  formulation  and  occurrence, 
it  also  proposes  a  solution.  Threats  targeting  high  valued  assets  such  as  a  carrier  would 
be  ranked  above  those  targeting  a  destroyer.  Similarly,  a  threat  arriving  earlier  would  be 
ranked  above  one  arriving  later.  But  what  makes  this  problem  not  straightforward  is  such 
concerns  as:  How  does  one  rank  a  threat  arriving  at  a  destroyer  earlier  than  one  arriving 
later  at  a  carrier?  The  ranking  scheme  proposed  here  is  based  on  the  “Loss  Threshold  Time” 
(LTT),  which  is  defined  as  the  time  at  which  the  expected  destructive  potential  of  a  threat 
exceeds  a  specified  value.  Then,  threats  with  earlier  LTTs  are  ranked  higher.  An  algorithm 
is  proposed  that  estimates  the  probability  of  each  defended  ship  becoming  the  target  of  a 
given  threat.  The  result  is  a  matrix  of  “Objective  Probabilities”,  which  is  used  to  calculate 
the  Expected  Loss”  (EL)  function  of  each  threat.  The  objective  probabilities  are  higher  for 
ships  that  can  be  reached  earlier  and  more  easily.  The  combination  of  these  two  concepts 
represents  a  conceptual  shift  in  how  warfare  planning  is  considered  in  algorithm  design  and 
implementation.  This  work  represents  a  significant  conceptual  breakthrough  that  may  have 
practical  consequences  in  future  combat  systems. 

This  work  was  supported  by  the  Office  of  Naval  Research  (Code  311). 

Applied  by:  f 

M&>L 

KENNETH  C.  BAILE,  Head 

Combat  Systems  Department 
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CHAPTER  1 
INTRODUCTION 


The  purpose  of  this  study  is  to  determine  for  each  individual  threat  in  an  AAW  or 
TBM  engagement  the  likelihood  that  the  combat  system  will  engage  a  defended  area  or  a 
particular  ship  in  a  task  force.  The  threats  include  both  Anti-Ship  Cruise  Missiles  (ASCM) 
and  Tactical  Ballistic  Missiles  (TBM),  Their  intended  targets  may  include  ships  at  sea  and/or 
extended  shore  areas  and/or  specified  shore  positions.  Any  of  these  potential  targets  are  to 
be  defended. 

Matters  of  engageability  or  fire  control  are  not  considered.  Each  threat  is  to  be  assigned 
a  value  indicating  both  the  potential  harm  it  is  capable  of  inflicting  and  the  urgency  by  which 
the  defenders  should  deal  with  it.  These  two  metrics,  of  “Potential  Harm”  and  “Urgency” 
are  quantized  so  that  they  are  available  to  the  decision-maker  that  prioritizes  and  engages 
the  threats.  Or,  these  two  measures  may  be  combined  into  a  single  value.  This  approach 
does  not  consider  the  decision  procedure  invoked  after  the  threat  evaluation. 

In  view  of  the  limited  scope  of  the  present  study  as  outlined,  the  inputs  required  in 
this  study  are  threat  state  information  and  the  positions  of  the  defended  points  and  areas. 
The  shooting  capacity  of  a  particular  defensive  platform  or  defense  strategy  is  not  needed. 
The  positions  of  the  ships  may  change  during  the  engagement.  An  assignment  of  value  is 
given  each  defended  point  and  area.  We  do  not  consider  complications  related  to  situational 
variation  in  value  (e.g.,  the  increased  value  of  a  surviving  carrier  after  loss  of  the  others). 

The  threat  information  ideally  consists  of  the  time  histories  of  all  threats.  The  past 
track  history  of  a  threat  undoubtedly  contains  information  of  use  in  its  identification  as 
well  as  prediction  of  its  potential  future  behavior.  However,  for  simplicity,  the  only  threat 
information  used  is  the  latest  threat-state  estimate  and  error  covariance  for  each  threat,  plus 
any  estimate  of  its  physical  parameters  that  may  be  available.  Usage  of  past  flight  history 
in  the  estimation  of  physical  parameters  is  assumed  to  have  been  carried  out  previously. 

The  metrics  Potential  Harm  and  Urgency  each  require  knowledge  or  a  hypothesis  re¬ 
garding  the  threat  goal.  Since  the  state  of  a  threat  cannot  be  known  with  certainty,  an 
“Objective  Probability”  is  estimated  for  each  of  its  possible  targets. 


1-1 


NSWCDD/TR-00/46 


Chapter  2  gives  an  overview  of  how  the  threats  are  evaluated.  Chapters  3,  4,  and  5  discuss 
the  treatment  of  ASCM  threats.  Chapters  6,  7,  and  8  discuss  the  TBM  threat.  Chapter  9 
discusses  the  method  by  which  the  Objective  Probabilities  and  other  evaluation  parameters 
are  estimated.  Chapter  10  discusses  the  organization  of  the  computations.  Chapter  11  gives 
a  summary  and  conclusions. 

Of  crucial  importance  is  the  estimation  of  the  Earliest  Time  of  Arrival  (ETA)  of  a  given 
threat  to  a  given  asset  (for  all  threat-asset  pairings).  Without  knowledge  of  the  intention  or 
capability  of  the  threat,  the  assumed  future  time  history  of  the  threat  position  is  hypothetical. 
The  assumptions  on  which  the  ETA  is  based  are  designed  to  result  in  a  simple  algorithm 
that  can  be  implemented  without  great  expenditure  of  computing  resources  and  still  provide 
a  usable  value  to  represent  a  rough  lower  bound  of  the  possible  ETA.  Thus,  the  objective  is 
to  represent  the  most  pessimistic  ( i.e .,  conservative)  situation.  The  ASCM  threat  is  assumed 
to  be  unpowered  and  flying  with  constant  speed.  Not  only  is  this  a  questionable  assumption, 
but  selection  of  the  speed  value  introduces  further  error. 

For  an  ASCM  this  value  is  set  equal  to  the  expected  speed  from  the  latest  state  estimate. 
For  a  small  to  moderate  threat  range  these  assumptions  may  not  be  too  bad.  But  when  the 
latest  observed  threat  altitude  is  high,  this  estimate  may  grossly  underestimate  the  average 
speed  and  correspondingly  overestimate  the  ETA.  This  suggests  a  future  upgrade  in  which 
a  “constant  energy”  rather  than  a  “constant  speed”  trajectory  is  used  in  some  cases  for  the 
estimation  of  ETA. 
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CHAPTER  2 

OVERVIEW  OF  THREAT  EVALUATION 


We  make  provision  for  treating  two  kinds  of  threats, 

ASCM 

TBM 

and  two  kinds  of  defended  positions, 

Defended  Points  ( e.g .,  ships) 

Defended  Shore  Areas  (e.g,,  cities) 

Treatment  of  ASCMs  differs  from  treatment  of  TBMs  in  the  following  respects.  Each  ASCM 
is  considered  potentially  able  to  target  any  of  the  Defended  Points.  These  are  given  “Ob¬ 
jective  Probabilities”  to  indicate  which  are  more  likely  to  be  targeted.  These  Objective 
Probabilities  therefore  depend  on  each  Threat-Objective  pairing. 

For  simplicity,  the  Objective  Probabilities  are  based  on  the  expected  state  of  the  threat, 
regardless  of  the  errors  in  the  state  estimation.  Maneuver  is  required  for  the  threat  to  reach 
any  of  the  Defended  Points.  The  Objective  Probability  is  high  for  points  that  can  be  reached 
quickly  (early  time  of  arrival)  and  easily  (low  drag  loss).  The  Objective  Probability  is  zero 
for  those  points  impossible  for  the  threat  to  reach  under  reasonable  assumptions  regarding 
its  flight  capabilities. 

In  contrast,  the  trajectory  of  a  TBM  is  considered  determined  by  its  characteristics 
(e.g.,  mass,  drag  coefficient,  etc.)  and  initial  state.  Maneuver  is  not  considered.  A  tra¬ 
jectory  calculated  by  simple  prediction  from  the  initial  expected  state  may  not  reach  any 
of  the  defended  assets.  But  these  are  possible  targets  because  the  initial  threat  state  and 
aerodynamic  characteristics  are  not  precisely  known.  Accordingly,  both  the  expected  ini¬ 
tial  threat  state  and  the  covariance  matrix  of  the  threat-state  estimation  errors,  as  well  as 
the  aerodynamic  uncertainties,  are  calculated.  These  determine  a  “footprint,”  which  is  the 
Probability  Density  Function  (PDF)  of  the  threat's  impact  point.  (The  reader  is  cautioned 
that  the  term  footprint  has  been  given  other  meanings.)  Defended  Points  lying  within  this 
footprint  are  assigned  Objective  Probabilities  related  to  the  footprint  PDF. 
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Extended  shore  areas  may  be  targets  of  TBMs,  but  not  of  ASCMs  unless  there  are 
Defended  Points  within  the  area.  The  Defended  Points  are  represented  by  their  positions, 
which  may  move  during  the  engagement,  and  by  an  “Actual  Value”  (Va),  and  a  “Perceived 
Value”  ( Vp ).  The  Perceived  Value  is  used  to  calculate  the  Objective  Probability,  while  the 
Actual  Value  determines  how  much  effort  is  to  be  expended  in  its  defense. 

A  Defended  Area  can  be  represented  geometrically  by  a  stationary  circle,  ellipse,  or 
a  polygon.  The  programs  and  examples  presented  here  use  ellipses.  A  TBM  generates 
an  Objective  Probability  determined  by  the  intersection  of  its  footprint  onto  the  Defended 
Area.  There  is  no  structure.  A  real  defended  area  such  as  a  city  might  have  variation  in 
value  throughout  its  area.  This  feature  can  be  simulated  by  placing  defended  points  within 
the  area. 

The  Objective  Probabilities  form  a  two-dimensional  matrix  whose  dimensions  equal  the 
number  of  threats  (active  at  a  given  time)  in  one  dimension,  and  the  number  of  Defended 
Points  and  Defended  Areas  in  the  other.  These  values,  along  with  the  Va  and  Vp  values, 
constitute  an  important  input  into  the  evaluation  of  the  threats. 

Engageability  of  the  threat  is  not  considered.  The  number  and  positions  of  the  shooting 
platforms  is  not  included  among  the  data  used.  No  PIPs  are  generated.  The  Objective 
Probabilities  are  estimated  without  regard  to  whether  or  not  the  threat  can  be  engaged,  or 
the  most  favorable  launch  times,  etc. 

2.1  Simplified  Trajectories 

For  each  (ASCM  -  Defended  Point)  pair  a  trajectory  producing  the  soonest  arrival  time 
is  determined.  The  initial  threat-state  error  covariance  is  ignored,  and  the  trajectory  is  based 
entirely  on  the  initial  expected  state.  The  initial  speed  is  maintained  constant  throughout 
the  flight.  The  trajectory  producing  the  earliest  arrival  time  is  determined.  This  will  usually 
lie  on  the  plane  determined  by  the  initial  threat  velocity  and  the  Defended  Point.  Initially 
there  is  a  maximum-g  turn  followed  by  a  CV  segment  to  the  Defended  Point.  Cases  where 
this  pattern  would  require  the  threat  to  fly  below  the  sea  surface  are  adjusted. 

The  initial  state  estimate  of  a  TBM  including  the  covariance  matrix  generates  a  suite  of 
trajectories,  one  for  each  possible  initial  state.  The  flight  maintains  constant  speed  with  no 
maneuver,  following  a  Keplerian  orbit  modified  for  liftless  drag  in  an  exponential  atmosphere. 

2.2  Determination  of  the  Objective  Probabilities 

At  any  given  time,  all  the  threats  that  are  active  at  that  time  are  considered  one  by  one. 
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The  ASCMs  are  treated  separately  from  the  TBMs. 

For  an  ASCM,  all  the  Defended  Points  are  considered  one  by  one.  For  each  Defended 
Point  the  drag  loss  along  the  simplified  trajectory  is  approximated  based  upon  available 
threat  information.  Incomplete  threat  information  is  assumed  to  satisfy  the  hypotheses  that 
it  minimizes  both  flight  time  and  drag  loss.  The  Objective  Probability  for  each  Defended 
Point  is  based  on  these  two  values.  After  all  Objective  Probabilities  have  been  approximated 
in  this  way,  they  are  normalized  so  that  their  sum  over  all  Defended  Points  is  unity. 

For  a  TBM,  the  first  step  is  to  estimate  its  footprint  in  order  to  determine  where  it  will 
land  on  the  Earth’s  surface.  The  expected  initial  threat  state  is  used  to  determine  a  single 
Kepler i an  trajectory  to  the  Earth  surface,  to  find  an  impact  position  and  arrival  time.  Next, 
the  atmospheric  effect  is  approximated  by  using  a  simplified  dynamic  model  of  a  mass  point 
falling  in  a  planar  orbit,  without  lift,  from  the  expected  initial  state  to  the  Earth’s  surface. 
The  resulting  flight  time  and  distance  provide  a  correction  to  the  previously-determined 
Keplerian  orbit  to  give  a  more  accurate  impact  position  and  flight  time. 

The  initial  state  covariance  is  assumed  to  be  Gaussian.  Even  with  this  idealization,  the 
complicated  nonlinear  flight  path  makes  it  difficult  to  get  a  closed-form  footprint  PDF.  In¬ 
stead,  the  impact  position  covariance  is  determined  by  assuming  small  errors  and  a  parabolic 
planar  trajectory  through  an  airless  medium.  This  impact  covariance  is  centered  on  the  ex¬ 
pected  impact  point  and  is  used  to  generate  a  two-dimensional  Gaussian  PDF  to  represent 
the  impact  footprint. 

Defended  Points  lying  within  the  impact  footprint  are  assigned  Objective  Probabilities 
related  to  the  Perceived  Value  and  the  PDF  density  at  that  position.  Defended  Areas 
intersecting  the  footprint  are  assigned  Objective  Probabilities  related  to  their  Perceived 
Values  and  the  integrated  footprint  probability  of  impact  in  the  intersection. 

2.3  Coordinate  Transformations 

Basically,  all  positions  are  referred  to  (rotating)  geographic  coordinates  (<f>,  L,  r).  Ve¬ 
locities,  on  the  other  hand,  are  represented  as  dimensional  components  ( u,v,w )  in  [L/T]  in 
the  directions  (<f>,L,r)  at  each  point.  (This  means  that  the  vector  represented  by  certain 
components  at  one  point  is  not  parallel  to  the  vector  represented  by  the  same  components 
at  another  point.)  It  is  useful  to  find  the  transformation  between  velocity  components  at 
two  different  points.  To  do  this,  it  is  convenient  to  define  a  Cartesian  coordinate  system 
(xhyu zi)  tangent  to  the  Earth  at  the  geographic  position  {<f>\ ,  Li , a),  with  x\  increasing 
Eastward,  yi  Northward,  and  zi  up.  This  will  be  referred  to  as  a  “Local  Cartesian  System”. 
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The  velocity  in  this  system  is  represented  by  the  components  in  the  directions 

A  velocity  vector  at  a  given  geographic  point  ( <j>,L,r )  is  represented  by  components 
(u,  v,w),  assumed  to  be  given  in  the  Local  Cartesian  System  referred  to  the  point  (<j),  L,  a) 
(directly  below).  It  is  required  to  relate  the  coordinates  of  points  and  the  components  of 
vectors  in  one  Local  Cartesian  System  to  those  in  another. 


Define  two  intermediate  coordinate  systems.  One  is  the  rotating  Earth-Centered  Carte¬ 
sian  System  (ECCS)  (xi,yi,zi),  where  y\  is  along  the  rotation  axis  pointing  toward  the 
North  Pole,  z\  is  in  the  direction  of  the  meridian  <f> i,  and  x\  ninety  degrees  East  of  this 
meridian.  The  transformation  between  (<f>,  L,  r)  to  (x\,y\,z\)  is 

x  i  =  r  cos  L  sin(^  —  <j>\ ) 

y\=rsmL  (2.1) 

z\  —  r  cos  L  cos  (<j>  —  <f>  i) 

The  second  is  a  rotating  Cartesian  Coordinate  system  (xi,yi,£i),  with  axes  parallel  to  the 
ECCS  but  centered  on  the  position  (<f>i,Li,a).  The  transformation  between  these  systems 
is 

A  /v 

X\  =  X\ 

V\  =  y\  —  a  sin  L\  (2.2) 

z\  =  z\  —  a  cos  L\ 

Then  the  system  (x\,y\,z\)  is  obtained  from  (x\,y\,z\)  by  a  rotation  about  the  x\  axis, 
according  to 

x\  =  x\  &i  =  x\ 

yi  =  yi  cos  L\  —  z\  sin  L\  y\  =  yj  cos  L\  +  z\  sin  L\  (2.3) 

z\  =  y\  sin  L\  +  z\  cos  L\  z\  —  —y\  sin  L\  +  z\  cos  L\ 

All  these  systems  are  non-inertial,  being  attached  to  the  rotating  Earth. 


Consider  now  a  second  Local  Cartesian  System,  and  corresponding  intermediate  systems 
defined  in  the  same  way  as  the  above,  with  all  the  subscripts  “1”  replaced  by  “2”.  Then 


transformation  of  the  components  of  velocity  referred  to  system  “1”  to  components  in  system 


“2”  is  done  by 


where 


U2\  d(x2,y2,Z2)  (  Ul\ 

v2  =  -XT - f  Vi 

w2  /  u{xi'>yiyZi)  / 


d(x2,y2,Z2) 

0(*i,yi,*i) 


(2.4) 


(2.5) 
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is  a  Jacobian  (matrix),  and  will  be  denoted  by  «/2l  •  This  matrix  can  be  found  by 

J21  =  J22  ^2g  ^ (2-6) 

where  the  suffix  g  stands  for  the  geographic  coordinates  ( cp,  £,  r).  Reversing  the  subscripts 
is  equivalent  to  taking  the  inverse  since  the  transformation  matrices  are  orthogonal.  It  is 
noted  that  and  Jj j  are  just  the  identity  matrices.  Since  J21  transforms  between  two 
Cartesian  systems,  its  value  is  independent  of  position.  This  can  be  seen  by  calculating  the 
product  of  the  inner  four  factors  of  (2.6)  obtaining 


Jp  ~  hlVglhx  — 


cos(^i  -(f) 2)  0  sin(^i  —  <f> 2) 

0  1  0 
sin(^2-^i)  0  cos(^i  -  <f>2) 


The  transformation  (2.6)  is  then 

(cos(^i  —  <f> 2)  — sinLi  sin(^i  —  ^2) 

sin  £2  sin(<^i  —  (f> 2)  cos  L\  cos  L2  +  sin  L\  sin  £2  cosf^i  —  fa) 

—cos  £2  sin(^i  —  (f>2 )  cos  £1  sin  £2  —  sin  L\  cos  £2  cos(^i  —  <f> 2) 

cos  £1  sin(^i  —  fa)  \  (^'®) 

sin  L\  cos  £2  —  cos  £1  sin  £2  cos(<fo  —  ^2)  1 
sin  L\  sin  £2  +  cos  L\  cos  £2  cos\fa  —  fa)  / 

As  an  example  of  the  importance  of  this  transformation,  consider  a  missile  at  30  degrees 
N  travelling  due  E  at  1000  m/s.  Relative  to  a  Local  Cartesian  System  based  on  a  point  one 
degree  further  East,  this  missile  has  a  Southward  velocity  of  8  m/s  and  an  upward  velocity 
of  15  m/s. 

2.4  Representation  of  Defended  Areas 

Defended  areas  are  represented  by  ellipses  with  specified  semimajor  axes  a,  semiminor 
axes  6,  centers  (20, 2/o),  and  directions  <f>a  of  their  lines  of  apsides  (with  <f>a  €  (— tt/2,  tt/2]). 
These  quantities  determine  the  eccentricities 


e  =  \/ 1  — 


and  the  positions  of  their  foci, 


Xfi  =  xq  —  aecosfa, 
Vfl  =2/0  —  aesin  <f>a, 

The  equation  for  the  ellipse  in  polar  form  is 


x/2  =  xo  +  ae  cos  <f>a 
Vf2  =  yo  +  ae  sin  <f>a 


a(l  -  e2) 

1  —  e  cos(0  —  <f>a) 


(2.10) 


(2.11) 
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where  (r,  9)  are  polar  coordinates  relative  to  the  point  Either  focus  could  have 

been  used.  For  (r,  9)  to  represent  polar  coordinates  relative  to  the  point  (z/2>y/2)>  the  first 
minus  sign  in  the  denominator  of  (2.11)  is  replaced  by  a  plus.  The  area  within  the  ellipse  is 

area  =  nab  (2.12) 

Each  defended  area  j  is  assigned  a  “perceived  value”  per  unit  area,  wj^  and  an  “actual 
value”  per  unit  area,  wja\ 


2-6 


NSWCDD/TR-00/46 


CHAPTER  3 

SIMPLIFIED  ASCM  TRAJECTORIES 


It  is  necessary  to  prioritize  engagements  against  threats  in  scenarios  with  a  large  num¬ 
ber  n  of  threats,  and  also  a  large  number  m  of  possible  targeted  assets.  The  amount  of 
computation  for  an  engagement  schedule  needs  to  be  kept  within  reasonable  bounds.  The 
objective  of  each  threat  is  not  known,  but  the  number  of  threat-objective  pairs  is  of  the 
order  of  nm.  Even  if  the  objective  were  known,  the  threat  trajectory  is  not.  If  all  possible 
trajectories  are  considered  the  amount  of  computing  can  easily  grow  unacceptably  large  for 
only  a  few  engagements.  For  this  reason,  the  “Simplified  Trajectories”  approach  presented 
in  this  section  is  needed. 

The  information  needed  to  prioritize  the  threats  includes  arrival  times  to  all  defended 
assets.  The  arrival  time  to  a  given  asset  for  a  given  threat  is  determined  by  a  single  simplified 
hypothetical  Earliest  Arrival  Time”  (ETA)  trajectory.  This  trajectory  is  based  on  the 
following  assumptions: 

•  A  maximum  allowed  lateral-acceleration  magnitude,  am,  is  specified.  This  value  can 
be  determined  by  intelligence  data  or  classification  based  on  imaging  or  past  threat 
behavior.  A  generic  value  is  used  if  there  is  complete  ignorance  of  the  threat  type.  The 
value  am  can  be  altered  at  any  time  additional  information  becomes  available. 

•  A  maximum  allowed  heading  error  (HE),  <j)m ,  is  specified.  The  HE  is  the  angle  between 
the  threat  velocity  and  the  asset  direction.  The  value  <j)m  attempts  to  represent  the 
limitation  on  the  threat’s  seeker  geometry.  Ideally,  we  would  be  concerned  with  the  look 
angle,  i. e.,  the  angle  between  the  threat  body  axis  and  the  asset  direction.  However, 
the  simple  analysis  proposed  here  treats  the  threat  as  a  moving  point,  and  consequently 
there  is  no  body  axis  to  which  to  refer  the  look  angle.  The  velocity  vector  serves  as  a 
reasonable  substitute,  and  thus  the  limitation  <f>m  is  placed  on  the  HE.  As  in  the  case  of 
am,  the  value  <f>m  can  be  determined  by  intelligence  data  or  classification  based  on  any 
available  information,  and  can  be  altered  during  the  course  of  the  battle. 

•  Constant  speed  is  assumed.  This  speed  is  determined  by  the  latest  threat-state  estimate. 
It  is  evident  that  these  assumptions  require  that  the  ETA  trajectory  lie  on  a  single  flight 
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plane  and  consist  of  two  segments,  “primary”  and  “secondary”.  The  primary  segment  is 
CTR  with  the  maximum-allowed  lateral  acceleration  am,  terminating  when  the  velocity  is 
directed  toward  the  targeted  asset.  The  secondary  segment  is  CV  to  reach  the  asset.  It  will 
be  shown  below  that  if  the  HE  constraint  is  initially  satisfied  it  will  remain  so. 

There  is  also  a  “Latest  Time  of  Arrival”  (LTA)  trajectory,  defined  as  a  single  CTR  arc. 
It  may  be  possible  to  artificially  reduce  the  threat  speed  by  applying  a  “maneuver  factor”, 
/m •  As  a  guide  for  selecting  /m,  note  that  a  planar  weaving  track  consisting  of  alternate 
CTR  turns  is  longer  than  the  corresponding  CV  track  by  the  factor  1  / fm  —  </>/sin (f>,  where 
<j>  is  the  maximum  heading  error.  Even  though  the  LTA  trajectory  is  not  used  at  this  time 
for  threat  evaluation,  it  is  included  in  our  discussion  both  for  general  interest  as  well  as  for 
providing  the  most  natural  limiting  case  for  delayed  arrivals. 

The  previous  paragraphs  define  the  two  limiting  trajectories  including  the  exceptional 
case  when  blind  application  of  the  trajectory  rule  drives  the  threat  trajectory  beneath  the- 
sea  surface,  as  may  happen  in  some  cases  with  a  diving  threat.  This  case  is  taken  up  in 
Section  3.3  and  in  Chapter  4. 

Let  fo,  rs,  and  fo  be  the  given  initial  threat  and  objective  positions  and  threat  velocity. 
The  vectors  rs  —  fo  and  ro  determine  the  plane  upon  which  the  trajectory  is  assumed  to  lie. 
On  this  flight  plane,  two  very  idealized  trajectories  are  constructed.  One  is  used  to  find  an 
approximate  ETA,  the  other  to  find  an  approximate  LTA.  The  former  (ETA)  is  an  important 
quantity  used  in  the  threat  prioritization.  The  latter  (LTA)  is  of  lesser  interest. 

3.1  The  ETA  Trajectory 

It  is  self-evident  that  the  trajectory  allowing  the  soonest  arrival  will  execute  a  maximum 
acceleration  magnitude  arc  until  headed  toward  the  objective  and  then  fly  CV  (constant 
velocity).  While  it  is  impossible  to  guess  the  threat’s  intention,  it  can  at  least  be  said  that 
the  arrival  time  determined  by  this  ETA  trajectory  provides  a  lower  bound  for  the  arrival 
time. 

Consider  first  the  maximum  g  arc.  Let  the  center  of  curvature  lie  at  the  point 

rc  =  fo  +  Afi o  +  B{f0  -  r,)  (3.1) 

where  A  and  B  are  constants  to  be  determined.  The  conditions  determining  these  constants 

|fo  -  re |2  =  b2m 

fo  •  (fo  -  fc)  =  0 
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These  conditions  determine  the  constants  in  (3.1)  to  be 


B  = 


±bmVt 


y/v? |ro  -  rs|2  -  [fo  •  (r,  -  r0)]2 


A=T-r(  r’-ro)B 


V? 


(3.3) 


If  the  ladical  in  (3.3)  is  too  small  for  practical  computation,  the  direction  of  tq  is  close  to 
the  direction  of  rs  —  ro,  (or  in  the  opposite  direction)  and  the  entire  trajectory  can  be  taken 
as  CV.  (A  test  should  be  performed  to  exclude  the  case  where  the  directions  are  opposite 
to  each  other,  which  normally  is  excluded  by  the  HE  constraint.)  The  two  signs  in  (3.3) 
determine  points  on  the  flight  plane  right  and  left  of  the  threat 5 s  directional  motion.  The 
point  lying  closer  to  the  objective  is  the  appropriate  arc  center,  and  is  determined  by  the 
condition 

(rc  -  ro)  •  (rs  -  ro)  >  0  (3.4) 

This  condition  selects  the  negative  sign  in  (3.3). 

Consider  the  point  rt  where  flight  transitions  from  the  arc  to  the  CV  segment.  Movement 
in  the  CTR  segment  can  be  represented  by  the  formula 


r(f)  =  ro  +  w  Vosincut-f  cv  Vo (1  -  coscut)  (3.5) 

(see  Groves  et  al,  1992  or  1994),  where  the  subscripts  “0”  refer  to  the  initial  state,  and 


u  = 


b 


[ro 


(3.6) 


being  the  turning  rate  and  b  the  turn  radius.  With  a  maximum  “g”  turn,  the  initial  accel¬ 
eration  is 


.  _  vM  -  f„) 

bl 


(3.7) 


with  bm  —  Vj* /  am  being  the  minimum  allowed  turn  radius.  If  t*  is  the  time  when  the  threat 
reaches  the  transition  point,  we  have 


The  condition 


requires  that 


ro  +  V  1bm  rosin  uiU  +  (1 


(rt  -  rc)  •  (ft  -  fs)  =  0 


C  cos  uitf  +  D  sin  uitt  =  bt 


(3.8) 

(3.9) 

13.101 


where 
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C  =  b2m-(fc-  f0)  •  (f,  -  ro),  D  =  Vf^mfo  •  (r,  -  r0) 


(3.11) 


If 


C2  +  D2  <  bl 


(3.12) 


there  is  no  solution.  Otherwise,  the  argument  of  the  trig  functions  in  (3.10)  is  given  by  either 
of 

ai  =  utt  =arccos  (b2Jy/ C2  +  £>2)  -  atan2(— D,  C) 

'  _ '  (3.13) 

a2  =  u>tt  —  2tt— arccos  (bm/V C2  +  ~D2  j  —  atan2(— Z),  C) 

These  values  are  placed  into  the  interval  [0,27r)  by  adding  or  subtracting  27r  if  necessary. 
The  smaller  value  is  then  used  to  determine  rt  from  (3.8).  At  this  point,  ft  is  tested  to 
insure  that  it  lie  above  the  sea  surface  (zt  >  0).  Otherwise  the  threat  is  unable  to  reach  its 
objective  by  the  simplified  trajectory  under  discussion  here,  and  the  possibility  of  a  different 
suitable  trajectory  is  considered  in  Section  3.3  and  Chapter  4. 

It  can  be  demonstrated  that  if  the  HE  constraint  is  satisfied  initially,  it  remains  satisfied. 
Consider 

f  •  (rs  -  f) 


COS  (f>  = 


(3.14) 


(3.15) 


V,\r,  -  v 

In  a  special  coordinate  system  in  which  ro  is  translated  to  the  origin,  rotated  so  that  rc  = 
(—bm,0)  and  the  linear  dimensions  scaled  by  defining 

xs  —  £/bn 
Vs  —  V/bn 

and  representing  the  (x,  y)  coordinates  of  a  point  on  the  CTR  arc  in  terms  of  the  polar  angle, 

x  =  cos  6 
y  =  sin  0 

(3.14)  becomes 

cos  0  = 


(3.16) 


— £  sin  6  +  r)  cos  6 


y/^2  +  V2  +  1  —  2£  cos  0  +  t)  sin  6 
The  initial  rate  of  change  of  (3.17)  with  respect  to  6  is 

±,r„M  -a?  ±  i2)  -  e 

d<>  \/[«  -  l)2  +  >)2P 

Since  for  a  valid  trajectory  the  point  (xs,ys)  must  lie  outside  the  circle 

(^s  +  bm)2  +  y%  —  b?„ 


(3.17) 


(3.18) 
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72>l-(£+l)2 


(3.19) 


must  hold.  Then  the  expression  (3.17)  is  always  positive  for  negative  £  showing  that  <f>  is 
always  increasing.  The  situation  at  positive  £  is  irrelevant  because  in  that  case  the  CTR  arc 
would  turn  in  the  opposite  sense. 

3.2  The  LTA  Trajectory 

It  is  difficult  to  find  an  objective  way  to  estimate  how  late  the  threat  could  delay  its 
arrival.  A  reasonable  assumption  is  that  he  will  not  want  to  delay  more  than  required  to 
execute  a  few  evasive  maneuvers.  In  any  case,  estimation  of  LTA  is  not  as  important  as  the 
ETA,  since  it  is  not  used  for  threat  evaluation. 

Define  a  rough  LTA  trajectory  as  consisting  of  a  single  circular  arc  on  the  same  flight 
plane  as  was  used  for  ETA.  This  will  provide  an  arrival  time  somewhat  later  than  the  ETA, 
but  perhaps  not  late  enough.  Further  massaging  may  be  done  by  subsequently  applying  a 
“Maneuver  Factor”. 

Again  representing  rc  by  (3.1),  determination  of  the  center  of  curvature  is  somewhat 
different  than  for  ETA  because  the  radius  of  curvature  b  is  not  known  a  priori.  The  conditions 
to  determine  rc  are 

\r,  ~  rc |2  =  |r0  -  rc|2  =  b2 

u  ~  -.s  n  (3-20) 

r o  •  (ro  —  rc)  =  0 

requiring  that  the  arc  radius  b  be  the  same  at  the  initial  and  final  threat  positions.  These 
conditions  lead  to  the  following  linear  equations  in  A  and  B. 


V%  A  +  ro  ■  (ro  —  rs)B  =  0 
2ro  •  (r0  -  rs)A  +  2lr0  -  rs\2B  =  — jro  -  rs| 


(3.21) 


After  B  and  A  are  determined  by  (3.21),  the  arc  radius  b  is  given  by  (3.14),  and  the  arc 
length  0S  is  given  by 

cos  6S  =  (r0  -  rc)  •  (fs  -  rc)/b2  (3.22) 

A  test  is  now  made  to  insure  that  the  LTA  trajectory  thus  obtained  lie  entirely  above  the 
sea  surface.  The  direction  of  the  derivative  (3.13)  with  6S  and  rs  replacing  0t  and  rt  will  be 
downward  in  a  dry  trajectory.  If  this  is  not  the  case  the  LTA  trajectory  can  be  modified 
to  keep  it  out  of  the  water  (see  Section  3.3).  Examples  of  the  ETA  and  the  LTA  trajectories 
are  shown  in  Figure  3.1 
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3.3  Keeping  the  ASCM  Above  the  Sea  Surface 

If  the  ETA  trajectory  as  described  above  indicates  flight  below  the  sea  surface,  an  al¬ 
ternative  “dry”  trajectory  is  sometimes  available  (see  Chapter  4).  If  the  ETA  trajectory 
is  entirely  above  the  sea  surface  but  the  LTA  trajectory  is  not,  an  alternative  “dry”  LTA 
trajectory  is  possible.  Whether  or  not  the  original  LTA  trajectory  crosses  the  sea  surface 
is  determined  by  examining  the  sign  of  the  vertical  velocity  component  upon  arrival  at  the 
target.  If  zs  <  0,  the  final  velocity  is  downward  and  the  original  LTA  trajectory  is  ok. 
Otherwise,  a  “dry”  LTA  trajectory  is  found  as  described  below. 


The  modified  trajectory  will  consist  of  an  initial  CTR  segment  followed  by  a  terminal 
sea-skimming  CV  segment.  The  point  of  transition  between  these  segments  is  denoted  by 
ft.  Both  segments,  and  the  transition  point,  lie  on  a  flight  plane  containing  the  points  ro 
and  rs,  and  oriented  so  that  it  contain  the  initial  velocity  vector  fo-  First,  determine  the 
hypothetical  “splash”  point  r;,  where  the  threat  would  hit  the  water  if  it  continued  in  CV 
flight  from  its  initial  state.  This  splash  point  is  given  by 


-  zoro 

n  =  r0 - — 

20 

(see  Figure  3.2).  It  is  evident  that  this  splash  point  also  lies  on  the  flight  plane.  From 
symmetry,  it  is  seen  that  the  distance  from  the  splash  point  to  the  transition  point  is  equal 
to  the  distance  from  the  splash  point  to  the  initial  position.  Call  this  distance  a.  Then  the 
transition  point  is  determined  by 


a(fs  -  fi) 

|rs  -  fi\ 


,  Vtz  o 

where  a  =  .  , 

M 


The  Pythagoras  relation  in  the  triangle  OTS  is 


(3.24) 


a2  +  62 


(3.25) 


determining  the  radius  of  curvature  b.  This  value  should  be  no  smaller  than  the  allowed  bm 
in  view  of  the  fact  that  the  ETA  trajectory  did  not  cross  the  sea  surface.  The  arc  length  of 
the  CTR  segment  is  then  given  by 


0  =  2arctan(a/6) 

and  the  length  of  the  CV  segment  is 


(3.26) 


(3.27) 
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SEA  SURFACE 


•  •  .F,lglfe  Geometry  of  the  ETA  and  LTA  Trajectories.  The  flight  plane  contains  the 
initial  threat  position  “0”  and  velocity  “V”,  and  defended  point  position  “S”.  The  two 
centers  of  curvature  are  indicated  by  “C”. 


c 

o 


o 


SEA  SURFACE 


Figure  3.2.  Geometry  of  the  “Dry”  LTA  Trajectory.  The  impact  point  “I”  is  where  the 
threat  would  hit  the  sea  surface  if  it  flew  an  initial  CV  track. 
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CHAPTER  4 

SEASKIMMING  TRAJECTORIES 

4.1  The  “Dual-Arc”  Conjecture 

This  conjecture  is  related  to  finding  an  alternative  simplified  ASCM  trajectory  when 
the  ETA  trajectory  determined  according  to  the  process  described  in  Chapter  3  is  partially 
under  water.  There  are  not  many  cases  where  this  will  occur  in  practice.  But  when  it  does, 
as  in  some  cases  of  diving  threats,  there  will  be  a  few  instances  where  a  valid  trajectory  can 
be  kept  out  of  the  water.  Although  this  situation  may  occur  rarely,  it  is  safer  to  account  for 
this  possibility  than  to  ignore  cases  where  an  asset  could  be  in  jeopardy. 

First,  a  definition: 

A  Valid  Trajectory  for  a  given  expected  initial  state,  targeted  point,  maximum  allowed 
heading  error  and  lateral-acceleration  magnitude  (ro,  ro,  rs,  (j)m,  and  a m)  is  one  consisting 
of  a  number  of  CV  and  CTR  segments  terminating  at  the  targeted  point,  with  continuity  of 
position  and  velocity  at  the  transition  points  between  segments,  where  the  speed  is  constant, 
and  z  >  0  everywhere. 

The  Conjecture:  If  a  valid  trajectory  exists,  there  will  also  exist  a  valid  trajectory  con¬ 
sisting  of  two  CTR  arcs,  the  final  arc  lying  along  the  sea  surface.  (Here  it  is  understood  that 
CV  is  a  special  case  of  CTR.)  We  have  been  unable  to  prove  this  and  invite  any  interested 
reader  to  communicate  his  results  to  us  if  he  is  able  to  find  a  proof  or  a  counter  example. 

Consider  a  flight  plane  of  a  primary  arc,  which  contains  the  initial  velocity  vector  0O  but 
otherwise  is  unconstrained,  allowing  a  one-parameter  family  of  such  planes.  A  Cartesian 
coordinate  system  (£,  r/,£)  is  attached  to  this  flight  plane  with  its  origin  at  the  initial  threat 
position  r0 .  Orthogonal  unit  vectors  (rgVj  ^ ,  u,  v)  form  a  right-handed  system  and  point 
in  the  direction  of  increasing  (£,?/,£)•  The  flight  plane  contains  the  vectors  F0  and  0,  is 
perpendicular  to  u ,  and  is  characterized  by  the  equation  rj  =  0. 

One  parameter,  A,  determines  the  orientation  of  the  flight  plane  according  to  the  relation 

— * 

u  =  a  cos  A  —  6  sin  A 

„  -  (4.1) 

v  =  a  sin  A  +  /?  cos  A 

where  (r0V~}  ,a,8)  is  another  right-handed  set  of  orthogonal  unit  vectors  that  obey  the 
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conditions 

a  •  ro  =  0,  /3  •  fo  =  0,  a  •  /?  =  0,  ft  ■  k  =  0,  a  -  k  >  0  (4.2) 

According  to  these  definitions,  the  unit  vectors  ci  and  $  are  related  to  the  Cartesian  unit 
vectors  ( i,j,k )  by 

—  (  nrr\0  _L  n  ir\  n  ^  _1_  1/2  V  _ .  0  /a')  Ta  O 

(4.3) 


a  = 


-(Aq?  +  yoj)zo  +  V£0k  ypT-  spy 

VSVfo  ’  1  Vh0 


where 


v^o  =  yjxl  +  yl 

A  point  on  the  flight  plane  having  coordinates  (£,  Q  is  located  at 

r  =  ro  +  roV^-1£  +  (,v 

In  other  words, 

f  =  {r  -  ro)  ■  roVf1 

C  =  (f  -  ro)  •  v 

The  Cartesian  components  of  (4.5)  are 

x  =  xo  +  4-  [-io^Osin(A)  +  yoVfCos(A)]Vrt_1Vrho1C 

y  —  yo  +  yoVt~l€  +  [-yoi*)sin(A)  -  ioVtcos(A)]V^_1Vr/lo1C 
z  =  zo  +  z0Vf^  +  V/l0sin(A)V^-1C 

Constant-speed  flight  in  a  circular  arc  can  be  expressed  in  the  form 

f(t)  =  fo  +  fbw'_1sin(u;^)  +  fow-2[  1  —  cos(u;t)] 

(see  Groves  et  al,  1992  or  1994)  where 

u;  =  |f|/|f| 

is  constant.  If  the  arc  (4.8)  is  tangent  at  the  sea  surface,  both  z  and  z  vanish  there, 
have 

z(tt )  =  zq  -f  ipu>-1sin(w£f)  +  3ow_2[l  —  cos(o>t*)]  =  0 
z(tt)  —  zocos(utt)  +  5o^-1sin(wti)  =  0 

at  the  instant  t  =  t*  the  threat  arrives  at  the  transition  point  ft  on  the  sea  surface, 
linear  combinations  of  (4.10),  it  is  seen  that 

/  .  \  .  /  j  \  “«DW 

cos  (bjit)  =  - — ; - 9,  sin  (u)tt)  =  — - 9 

ZQ  +  ZQU 1  ZQ  +  ZQUJ 1 


(4.4) 


(4.5) 


(4.6) 


(4.7) 


(4.8) 

(4.9) 
and  we 

(4.10) 
Taking 

(4.11) 
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Summing  the  squares  of  (4.11)  it  is  seen  that 

Zqu2  =  zl  -  2z0z0  (4.12) 

is  a  condition  that  an  arc  on  the  flight  plane  meet  the  sea  surface  at  a  point  of  tangency. 

Consider  the  problem  of  determining  r0  from  a  given  r0,  f0,  and  A.  The  conditions  (4.9), 
(4.12),  and  the  relation  rb  ■  u  =  0  give  the  equations 

%oxo  +  yoyo  +  hh  =  0 

UxX  o  4-  Uyjjo  4-  uzzq  =  0  (4-13) 

4(4 +4  +  4)  +  zqzo  =  zlv? 

Eliminating  xq  and  yo  from  (4.13)  gives  a  quadratic  equation  in  zq  that  may  or  may  not 
have  real  solutions. 


If  the  initial  lateral  acceleration  magnitude  |f0|  (and  hence  the  radius  of  curvature  b  = 
Vt  /|ro|)  is  given  instead  of  A,  then  u  is  unknown  and  instead  of  (4.13)  the  initial  acceleration 
is  determined  by  solving  the  equations 

xqxo  +  yoyo  +  zq'zq  =  0 

2 vt  zozo  =  4vt2  ~  ^ol^oi2  (4.14) 

4  +  4  +  4  =  In)  |2 


Consider  (4.14)  in  the  space  of  the  acceleration  components.  The  three  equations  rep¬ 
resent,  respectively,  a  plane  through  the  origin,  a  plane  parallel  to  the  (xo,yo)  plane,  and  a 
sphere.  The  intersection  of  the  plane  through  the  origin  with  the  sphere  occurs  at  values  of 
zq  contained  in  the  interval 


2o|  < 


Yh  o|n)| 
Vt 


(4.15) 


If  the  value  zq  determined  by  the  second  equation  of  (4.14)  lies  in  this  interval,  then  a  real 
solution  exists.  This  condition  requires  that 


Vt-Vh0  Vt  +  Vh0 

OJ  — - 


Z0  ZQ 

The  equation  in  xq  obtained  by  eliminating  yQ  from  (4.14)  is 

Vho4  +  2xzz0xq  +  ZqO/2  +  z2)  -  |f0| 24  =  0 


(4.16) 


(4.17) 


This  equation  has  real  roots  (for  x0 )  if  (4.16)  is  satisfied.  Then  the  flight  plane  (A)  is  obtained 
by  getting 

v  —  r0/|f0| 

cos  X  —  (5  •  v  (4.18) 

sin  A  =  a  •  v 
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Consider  all  the  possible  flight  paths  that  could  serve  as  the  primary  arc  without  consider¬ 
ation,  for  the  moment,  of  the  position  of  the  targeted  point.  Figure  4.1  is  a  three-dimensional 
plot  on  an  arbitrary  scale  of  such  an  example  in  which  a  one-dimensional  family  of  candidate 
trajectories  is  shown  leaving  the  initial  threat  position.  The  trajectories  displayed  corre¬ 
spond  to  values  of  u  at  equal  increments  between  the  limits  (4.16).  Each  arc  is  tangent  to 
the  sea  surface  at  a  point  indicated  by  a  dot.  These  points  of  tangency  lie  on  a  circle,  as  will 
now  be  demonstrated. 

Let  us  introduce  a  “Special  Coordinate  System”  (SCS)  on  the  sea  surface,  with  coordi¬ 
nates  (£,7/),  so  that  the  SCS  origin  is  under  the  initial  threat  position,  with  the  £  axis  in 
the  direction  of  the  initial  horizontal  velocity.  With  respect  to  the  heretofore  used  “Basic 
Coordinate  System”  (BCS)  on  the  sea  surface  with  coordinates  (x,y),  the  two  systems  are 
related  by 

£=  (x  —  a:o)cos  A  +  (t/ —  yo)sin  x  ~  x0  —  £cos  A  —  7/ sin  A 

7/  =  -(x  -  xo)sin  A  +  (y  -  yo)cos  A,  y  -  yo  =  £  sin  A  +  7;  cos  A 

where  A  is  the  azimuth  of  the  initial  horizontal  threat  velocity  in  the  BCS, 

A  =  atan2(r»o,7to) 


(4.19) 


(4.20) 


Equations  (4.10),  (4.11),  and  (4.12)  determine  the  coordinates  of  the  point  of  tangency  in 
the  SCS,  giving 


-to  +  £0*0  VQZQ 

..  9  j  7 Jt  ..  .  9 

ZQ  +  ZQU2  ZQ  +  ZQU)Z 


(4.21) 


It  is  convenient  to  use  uj  as  the  parameter  that  selects  one  out  of  the  family  of  arcs.  The 
initial  acceleration  is  determined  by  (4.9),  (4.12),  and  the  orthogonality  of  the  acceleration 
and  velocity  vectors.  The  result  is 


£0  =  to(*o  -  ZqJ2)/2zo^o 

7)0  =  ±Vt\[u2  -  (zq  -  zgw2)2/4zo£o 

zo  =  (zq  ~  ZqU2)/2z0 


(4.22) 


Substituting  these  values  into  (4.21)  we  obtain  the  equation  of  the  curve  passing  thru  the 
points  of  tangency  on  the  sea  surface  in  terms  of  the  parameter  u.  Eliminating  u>  by  first 
expressing  it  in  terms  of  £*, 


j2  _  ifl(2zo£o  +  zqZq  +  £qzq&) 
Zq{zqZ0  -  £o&) 


(4.23) 
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after  some  rearrangement  gives 


Vt  + 


(4.24) 


demonstrating  the  circular  nature  of  the  locus  of  tangency  points,  and  determining  the  center 
and  radius.  Note  that  the  center  of  this  tangency  circle,  according  to  (4.24),  is  just  the  point 
on  the  sea  surface  wheie  the  threat  would  impact  if  it  continued  in  a  CV  trajectory  from  its 
initial  condition. 


Just  as  amusing  is  the  fact  that  the  velocity  at  the  points  of  tangency  is  directed  outward 
from  the  center  of  this  “tangency”  circle.  This  velocity  is 

ft  =  rocos  utt  4-  roa/-1sin  cut*  (4.25) 

according  to  (4,8).  Combining  this  with  the  conditions  (4.11)  gives 

%£o  ~  ...  20*70 


&  = 


«  ,  9  5  —  7.  ~  a 

20  +  Z()UJ  ZQ  +  ZQUJ* 

Using  the  relations  (4.22)  and  (4.23),  it  can  be  shown  that 


(4.26) 


K  =  <4-27> 

where  p  is  the  vector  from  the  center  of  the  tangency  circle  to  the  point  of  tangency,  given 

by 

zoro 

P  =  r,  +  —  (4.28) 

These  vectors  are  shown  on  an  arbitrary  scale  on  the  plane  of  the  sea  surface  in  Figure  4.2  for 
the  arcs  of  Figure  4.1,  where  the  initial  threat  position  projected  down  onto  the  sea  surface 
is  indicated.  The  points  of  tangency  are  indicated  by  the  dots. 


Returning  to  the  BCS,  a  targeted  point  rs  is  added  to  the  example  of  Figure  4.1.  To 
find  the  secondary  segment  of  the  dual-arc  trajectory,  use  the  sea  surface  as  the  complex 
plane,  since  the  segment  sought  lies  on  this  plane.  To  find  the  arc  center  for  the  segment 
continuing  from  any  transition  point  ft,  note  that  its  complex  value  is 


Wc  =  Wt  +  bei&+*M  (4.29) 

where  b  is  the  (signed)  radius  of  curvature  and  Wt  is  the  position  of  the  transition  point, 
and  (f>  is  the  argument  (direction)  of  the  final  velocity  of  the  primary  segment.  The  radius  b 
is  to  be  determined  by  the  condition 


wc-wtf  =  \wc-ws\2 
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This  equation  gives  the  value 

t  W8W*  +  WtW*  -  WtW*  -  WsWt* 

~  {W*  -  Wt*)e^+7r/2)  +  (Ws  -  Wt)e-'(<t>+*/V 
_  VtVa-rtf 
2| ft  x  (rs  -  ft)| 

Since  both  factors  in  the  cross  product  of  (4.31)  lie  in  the  (x,  y )  plane,  the  product  is  in  the 
±k  direction,  where  k  is  in  the  upward  direction.  If  the  product  is  in  the  +k  direction  the 
turn  is  toward  the  left,  and  vice  versa.  The  length  of  the  arc  is  given  by 

6C  =  arccos[b-2(Ff  -  rc)  •  ( rs  -  fc)}  (4-32) 

Points  along  this  secondary  arc  are  conveniently  computed  using  equation  (4.8). 


Consider  the  limiting  simplified  trajectory,  either  ETA  or  LTA,  of  Chapter  3,  whose 
transition  point  is  on  the  sea  surface.  This  trajectory  will  lie  on  a  single  flight  plane,  while 
the  dual-arc  trajectories  generally  lie  on  two  flight  planes.  To  determine  the  single-plane 
limiting  trajectory,  find  the  “impact  point”  r*  on  the  sea  surface  if  the  threat  were  to  fly  CV 
from  its  initial  position.  This  point  is  given  by 


n  =  ro  - 


(4.33) 


Referring  to  Figure  3.2,  it  is  evident  from  symmetry  that  the  transition  point  ?t  is  given  by 


h  =  n  + 


(4.34) 


Figure  4.3  shows  a  plan  view  of  the  single  family  of  dual-arc  trajectories  as  solid  curves 
on  an  arbitrary  scale,  plus  the  limiting  simplified  trajectory  shown  by  the  dashed  curve.  The 
trajectories  have  been  projected  vertically  onto  the  sea  surface.  The  dots  indicate  the  points 
on  the  dual-arc  trajectories  where  the  primary  segment  is  tangent  at  the  sea  surface.  The 
portions  of  all  the  trajectories  lying  outside  these  tangent  points  are  on  the  sea  surface,  while 
the  portions  inside  these  points  are  above  the  sea  surface.  The  secondary  arc  is  shown  only 
if  the  heading  error  at  the  transition  point  does  not  exceed  7t/2. 

If  a  dual-arc  trajectory  is  needed,  the  dashed  limiting  trajectory  must  require  a  lateral 
acceleration  in  excess  of  the  allowed  maximum  value.  It  is  apparent  in  Figure  4.3  and  similar 
cases,  that  there  is  a  dual-arc  trajectory  (one  of  the  solid  curves)  in  which  the  maximum 
lateral  acceleration  is  less  than  that  in  the  single-segment  (dashed)  curve.  Thus,  there  exist 
valid  dual-arc  trajectories  in  some  cases  for  which  the  usual  simplified  ASCM  trajectory 
would  require  the  threat  to  fly  under  water. 


4-6 


NSWCDD/TR-00/46 


Figure  4.4  shows  in  perspective  on  an  arbitrary  scale  a  family  of  possible  primary  arcs 
for  an  initial  threat  velocity  that  is  somewhat  downward.  The  dots  at  the  ends  of  these  arcs 
lie  on  the  tangency  circle  on  the  sea  surface. 

Examples  of  other  possible  dual  arc  trajectories  are  shown  on  an  arbitrary  scale  in  Figures 
4.5,  4.6,  4.7,  and  4.8. 

4.2  Dry  Trajectories 

Consider  the  Tangency  Circle  Coordinate  System”  (T CCS)  in  which  the  origin  of  the 
coordinates  ( X ,  Y)  is  at  the  center  of  the  circle,  and  the  targeted  ship  is  at  (X,  Y)  =  (|rs  - 
rc|,0).  The  Y  axis  of  the  TCCS  intersects  the  f  axis  of  the  SCS  at  an  angle  a  given  by 

<7  =  atan2[(£s  -  fc) ,  -r;s]  (4,35) 

where  rc  is  the  center  of  the  tangency  circle  and  rs  is  the  position  of  the  targeted  asset.  The 
transformation  between  the  SCS  and  the  TCCS  is 

£  —  fc  =  Y  cos  er  +  .X"  sin  cr,  X  —  —7]  cos  a  +  (£  —  £c)sin  <7 

rj  =  Y  sina  -  X  cos  a,  Y  =  7;  sin  a  +  (f  -  fc)cos  a  {^3$) 

The  relations  between  the  three  coordinate  systems  are  shown  in  Figure  4.9.  Starting  with 
the  original  BCS  represented  by  the  (x,y)  axes,  the  SCS  is  determined  by  the  position  of 
the  threat  at  point  T  with  horizontal  velocity  in  the  direction  of  C.  The  £  axis  is  along  the 
line  TC,  which  makes  the  angle  A  with  the  x  axis,  with  the  point  C  being  determined  by  the 
offset  of  the  center  of  the  tangency  circle  indicated  in  (4.24),  The  SCS  origin  is  at  the  point 
T.  The  TCCS  then  has  its  origin  at  the  tangency  circle  center,  C,  with  the  X  axis  toward 
the  targeted  asset  position  S.  The  Y  axis  intersects  the  £  axis  at  an  angle  a.  The  locus  of  the 
secondary  arc  centers  is  the  line  parallel  to  the  Y  axis,  as  will  be  shown  below,  intersecting 
the  X  axis  at  the  point  P.  The  location  of  P  is  given  by  (4.42). 

We  can  assume  that  the  targeted  asset  lies  outside  this  tangency  circle.  Otherwise,  either 
a  valid  simplified  trajectory  would  reach  the  asset,  or  no  valid  dual-arc  trajectory  would  reach 
it.  In  either  case,  there  is  no  need  for  a  dual-arc  trajectory. 

Position  on  the  tangency  circle,  ra,  is  specified  by  the  angle  (3  measured  counterclockwise 
from  the  direction  of  the  targeted  asset,  according  to  the  relation 

Fa  (/?)  =  rc  +  bc(i  cos  /5  +  J  sin  (3 )  (4.37) 

where  here  i  and  J  are  unit  vectors  in  the  X  and  Y  directions,  and  bc  is  the  radius  of 
the  tangency  circle.  Recall  that  the  threat  velocity  at  a  point  on  the  tangency  circle  is 
perpendicular  to  the  circle. 
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The  secondary  arc  begins  at  ra{ft)  and  ends  at  r,  in  a  circular  arc.  Let  rac(ft)  be  the 
center  of  this  arc.  The  relations 


i?»  -  =  w)  -  f«m = h(P) 

-  ?„(/?)]  •  |r„(/3)  -  Fc|  =  0 


(4.38) 


determine  rac(ft),  where  62  is  the  (unspecified)  arc  radius.  Let 


rac  —  fc  +  Ara  4-  Brs 


(4.39) 


where  A  and  B  are  to  be  determined  by  (4.38). 

In  terms  of  the  TCCS  coordinates  X  and  V,  (4.38)  becomes 

2 bc{bc  -  Xscosft)A  +  2Xft(bccos  ft  -  XS)B  =  b2c-X2 

Ab^.  "f  B X$ cos  ft  =  be 


(4.40) 


A  = 


giving 

2bgXs  -  (b2  +  Xl)cosft 
2bcXss\r?ft 

The  locus  of  the  secondary  arc  centers  is  then 

b2  +  X f 


B  = 


l?c-2bcXscosft  +  X'> 
2X|sin2/5 


Xac  ~ 


2X, 


2bcXs-(b2c  +  Xi)cosft 
ac  2Xssin  ft 


(4.41) 


(4.42) 


according  to  (4.39).  It  is  seen  that  this  locus  is  a  line  because  Xac  is  independent  of  ft. 

The  geometry  of  this  secondary-arc  construction  is  shown  in  Figure  4.10  as  a  family  of 
possible  arcs.  The  dashed  circle  is  the  tangency  circle.  The  threat  position  is  generally  not 
at  the  center  of  this  circle.  The  dashed  line  is  the  locus  of  arc  centers. 


Consider  the  region  of  the  tangency  circle  permitted  by  “g”  constraint.  For  the  moment 
we  shall  ignore  the  constraint  on  Heading  Error  and  consider  only  how  the  limitation  on 
the  maximum  allowed  lateral  acceleration  magnitude  limits  the  possibilities.  Representing 
(4.16)  in  terms  of  the  radius  of  curvature,  b\,  of  the  primary  arc,  there  is  the  relation 


ba  <  b\  <  b\, 


(4.43) 


where 

U  _  zoVt  L  __  zoVt 

°  Vt  +  Vh o’  *  Vt-VM 


(4.44) 


If  bm  <  ba  the  simple  trajectory  will  generate  a  transition  point  r)  above  the  sea  surface, 
and  there  is  no  need  to  consider  the  dual-arc  trajectories.  If  65  <  bm  the  threat  cannot  turn 
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sharp  enough  to  avoid  the  sea  surface  regardless  of  how  he  orients  his  flight  plane,  and  thus 
no  dual-arc  trajectory  is  possible.  The  remaining  case  is  where 


ba  5;  bm  Y  bjj 


(4.45) 


There  always  exist  valid  trajectories  for  the  secondary  arc,  because  the  arc  from  the 
point  on  the  tangency  ciicle  nearest  to  rs  is  a  straight  line.  Moving  along  the  tangency  circle 
in  eithei  direction  the  associated  arcs  reach  a  minimum  radius  of  curvature,  then  again  get 
larger  without  limit.  The  arcs  with  minimum  radius  of  curvature  have  their  center  on  the 
line  between  rc  and  rs.  The  radii  of  these  arcs  is  given  by 


bo,m  — 


1^  ~  rc]2  -  b\ 
2|r«  -  rc | 


(4.46) 


This  is  the  minimum  radius  of  all  the  possible  secondary  arcs.  If  bm  <  b2m  the  secondary  arcs 
departing  from  any  point  on  the  tangency  circle  are  valid.  If  b2m  <  bm  some  are  not  valid. 
To  find  the  points  (X2m,Y2m)  on  the  tangency  circle  corresponding  to  the  arc  of  minimal 
radius,  the  conditions 


(^2 m  ~  X ac )2  +  Y2m  =  b\ 


2m 


X2 m  +  Y2m 


?  —  }? 


(4.47) 


determine  two  points  on  the  tangency  circle. 


The  situation  is  described  by  Figure  4.11,  which  again  depicts  the  plane  of  the  sea  surface. 
The  point  “T”  is  the  projection  of  the  threat  position,  “S”  is  the  targeted  ship,  and  “C”  is 
the  center  of  the  tangency  circle.  If  condition  (4.43)  holds,  there  is  a  segment  of  the  tangency 
circle  associated  with  valid  trajectories  in  the  primary  arc.  The  points  “fa”  and  “fa”  denote 
the  points  of  tangency  of  primary  arcs  having  the  limiting  radius  of  curvature  bm,  and  are 
determined  by  (4.21)  and  (4.22).  If  the  maximum  allowed  turning  rate,  =  Vtfbm,  satisfies 
(4.16),  then  the  two  final  primary-arc  positions  on  the  tangency  circle  are  determined.  The 
corresponding  values  fa  and  fa  divide  the  tangency  circle  into  a  permitted  and  a  prohibited 
arc.  The  permitted  arc  is  identified  by  the  fact  that  it  must  include  the  angle  fa  =  7t/2  -  a 
(see  Figure  4.9). 


The  dashed  line  perpendicular  to  CS  is  the  locus  of  arc  centers  for  possible  trajectories 
in  the  secondary  arc.  The  points  “Pi”  and  “P2”  are  where  the  secondary  arc  of  minimum 
radius  intersects  the  tangency  circle,  and  has  coordinates  (X2m,  Y2m )  as  determined  by  (4.47). 
There  are  two  solutions,  which  determine  the  two  points.  The  center  of  this  arc  is  the  point 
“P”.  if  b2m  >  bm,  then  all  the  secondary  arcs  satisfy  the  “g”  constraint.  If  &2m  <  bm  the 
secondary  arcs  arriving  at  the  tangency  circle  in  the  vicinity  of  the  points  “PI”  or  “P2” 
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require  a  turn  sharper  than  allowed.  In  this  case,  to  find  the  secondary  arc  having  the  radius 
of  curvature  bm,  consider  the  conditions 


(Xm  ~  Xacf  +  (Ym  -  Ycl)2  =  b2m 


Xl  +  Y2=b l 


(4.48) 


where  Yc\  is  the  (positive)  coordinate  of  the  point  Cl  in  Figure  4.11,  given  by  the  relation 


Y?i  +  (X,  -  Xmf  =  bi 


(4.49) 


If  there  is  no  real  solution  for  Yc\,  then  all  the  possible  secondary  arcs  satisfy  the  ‘g’  con¬ 
straint.  Otherwise,  there  are  two  solutions  to  (4.48),  (Xmi,Yml)  and  (Xm2,  Ym2),  with  Ym l 
and  Ym 2  having  opposite  signs.  Define  Ym\  as  the  negative  solution.  Then  the  azimuth 
angles  on  the  tangency  circle  where  these  secondary  arcs  impinge  are  given  by 


71  —  —  &tan2(V^ii,  Xm\) 

72  =  atan2(V^l2>  Am2) 


(4.50) 


with  both  values  being  positive.  Secondary  arcs  arriving  at  the  tangency  circle  within  the 
arc  “—72  to  —71”  or  within  “71  to  72”  have  radii  less  than  allowed,  and  accordingly  are  not 
valid.  The  remaining  points  along  the  tangency  circle  have  secondary  arcs  that  satisfy  the 
‘g’constraint. 


Consider  the  region  of  the  tangency  circle  permitted  by  HE  constraint.  Figure  4.12 
shows  the  tangency  circle  and  the  plane  of  the  sea  surface  including  the  targeted  ship  “S” , 
the  projection  of  the  threat  position  “T” ,  the  center  of  the  tangency  circle  “C” ,  and  the  line 
(dashed)  locus  of  the  secondary-arc  centers.  Consider  a  secondary  arc  whose  center  is  denoted 
by  “Cm” .  It  is  evident  that  the  HE  along  this  arc  is  greatest  where  this  trajectory  departs 
from  the  tangency  circle.  If  this  angle  is  equal  to  (f)m,  the  maximum  allowed  HE,  denote 
the  azimuth  angle  P  on  the  tangency  circle  where  these  limiting  secondary  arcs  intersect  the 
circle  by  — /3m  and  Pm. 


Denote  the  TCCS  coordinates  of  the  lowest  of  the  pair  of  critical  points  by  (Xm,Ym), 
given  by 

Xm  =  bccospm,  Ym =  bc  sin  Pm  (4.51) 


where  bc  is  the  radius  of  the  tangency  circle.  Note  that  the  angle  at  S  subtended  by  the 
points  C  and  —  pm  is  equal  to  (j)m  —  Pm-  Accordingly,  there  is  the  relation 


tan  ((/>m  —  Pm)  — 


be  sin  Pm 

Xs  be  COS  pm 


(4.52) 


or 


0m  =  Pm  +  atan2  (bc  sin  pm,  Xs  -  bc  cos  Pm) 
4-10 


(4.53) 


NSWCDD/TR-00/46 


Taking  the  tangent  of  both  sides  and  rearranging  gives  a  quadratic  equation  in  sin  /3m, 

Xg  sec2 (pms'm2 (3m  +  2Xsbctan  <f)m sin  pm  +  (b2c  -  X2)tan2(/)m  =  0  (4.54) 

The  discriminant  of  this  equation  is 

\J 4*f  ten ~  6i)tan20m  +  Af]  (4.55) 

and  is  always  real  for  a  targeted  ship  outside  the  tangency  circle  (i.e.,  Xs  >  bc).  It  is  also 
seen  that  the  largest  root  is  always  non-negative.  After  obtaining  sin  j3m  from  (4.54)  the 
appropriate  quadrant  of  /3m  is  found  using  (4.53),  giving  /3m  on  [0.7r). 

The  arc  from  the  point  labeled  /?m  to  S  is  obtained  by  symmetry.  The  permitted  region 
of  the  tangency  circle  under  the  HE  constraint  is  the  arc  between  —  /3m  and  f3m. 

Consider  now  the  HE  constraint  on  the  primary  arc.  The  geometry  is  complicated  in  the 
general  case,  and  therefore  the  following  simplifying  assumption  is  made.  This  assumption 
can  be  phrased  in  the  form  of  a  conjecture: 

If  both  the  initial  and  final  HE  of  the  primary  arc  satisfy  the  HE  constraint,  the  HE  at 

all  points  of  the  arc  satisfy  the  constraint. 

It  may  be  possible  to  show  that  this  conjecture  is  true,  although  we  have  been  unable 
to  do  so.  However,  if  there  are  some  unlikely  cases  where  the  HE  constraint  is  violated 
somewhere  along  the  primary  arc  even  though  the  constraint  is  satisfied  at  the  end  points, 
the  HE  there  will  only  slightly  exceed  the  allowed  limiting  value.  The  HE  at  the  end  of  the 
primary  arc  is  the  same  as  at  the  beginning  of  the  secondary  arc.  Therefore,  the  permitted 
region  of  the  tangency  circle  is  the  arc  /3  6  [/?i,  /%]• 

Figure  4.13  shows  how  the  geometry  of  Figure  4.11  determines  the  regions  of  validity.  A 
point  on  the  tangency  circle  determines  both  the  primary  and  the  secondary  arc.  Regions  of 
this  tangency  circle  are  mapped  into  the  angle  6 ,  measured  counterclockwise  around  the  circle 
from  the  point  nearest  to  “S”.  In  Figure  4.12,  fi  along  the  tangency  circle  has  the  values 
indicated  by  the  point  labels.  Figure  4.13  depicts  a  situation  with  the  most  complicated 
region  of  valid  trajectories,  namely,  three  disjoint  intervals.  Other  situations  will  correspond 
to  two  intervals  or  a  single  interval. 
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Figure  4.2.  Points  of  tangency  on  the  sea  surface  of  the  arcs  shown  in  Figure  4.1,  indicated 
by  the  dots.  The  line  segments  indicate  the  velocity  direction  at  these  points. 
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Figure  4.6.  Projection  onto  the  sea  surface  of  two-arc  trajectories  reaching  the  indicated 
targeted  point  (solid  curves).  The  initial  arcs  are  those  of  Figure  4.4  departing  from  threat 
located  above  the  origin.  The  dashed  curve  is  the  trajectory  lying  on  a  single  flight  plane, 
terminating  in  CV  segment.  •  ' 
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Figure  4.7.  Projection  onto  the  sea  surface  of  two-arc  trajectories  reaching  the  indicated 
targeted  point  (solid  curves).  The  initial  arcs  are  those  of  Figure  4.4  departing  from  threat 
located  above  the  origin.  The  dashed  curve  is  the  trajectory  lying  on  a  single  flight  plane, 
terminating  in  CV  segment. 
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Figure  4.9.  Relation  between  the  BC  (z,j/),  SC  and  TSC  {X,Y)  Systems. 
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Figure  4.12.  The  tangency  circle  and  limiting  secondary  arcs  determined  by  the  HE  con¬ 
straint. 
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Figure  4.13.  Regions  of  validity  along  the  6  axis. 


NSWCDD/TR-00/46 


CHAPTER  5 

REGIONS  OF  AS  CM  ACCESSIBILITY 


In  accordance  with  the  assumptions  defining  the  simplified  ASCM  trajectory,  the  threat 
is  not  allowed  to  reach  points  whose  heading  error  exceeds  (f>m ,  nor  requiring  a  turn  of  radius 
less  than  bm.  The  first  requirement  restricts  the  Region  of  Accessibility  (RA)  to  the  interior 
of  a  circular  cone  whose  vertex  is  at  the  current  expected  threat  position  and  whose  axis  is 
directed  along  the  expected  threat  velocity.  The  second  condition  excludes  from  the  RA  the 
interior  of  a  torus  to  be  described  below. 


It  should  be  emphasized  that  the  current  RA  relates  to  the  current  threat  state.  One 
might  think  that  the  RA  at  a  subsequent  time  would  be  contained  within  the  earlier  RA, 
but  this  is  not  necessarily  so.  At  any  time  the  RA  only  represents  the  region  accessible  to 
simplified  trajectories  determined  by  the  threat  state  at  that  time. 


Consider  the  conical  region.  A  point  rs  is  in  the  allowed  conical  region  of  the  RA,  as  in 
Figure  5.1,  if 

r  •  (r  —  f ) 

_  r|  —  cos  C5-1) 

where  r  =  (x,y,  z)  and  r  =  (u,  v,  w)  are  the  current  expected  position  and  velocity  of  the 
threat  and  V  =  |r|.  For  a  point  on  the  sea  surface,  rs  =  (xs,ys,Q),  the  condition  (5.1) 
imposes  the  constraint 


A£i  +  B£ \T)i  +  Ciji  +  D&  4-  Erji  +  F  >  0 

where 

m  =  ys  -y 

and 

A  =  u2  -  V2cos2(/>m  D  =  —2 uwz 

B  =  2  uv  E  =  -2  vwz 

C  =  v2  —  V2cos24>m  F  =  z2(w 2  —  V2  cos2  (pm) 

The  terms  linear  in  and  r/i  can  be  eliminated  by  the  substitution 

€i=a£  +  bri  +  e 
m  —  +  dr]  +  / 


(5.2) 

(5.3) 


(5.4) 


(5.5) 


(5.6) 


and  setting 


2Aae  4-  B(af  4-  ec)  4-  2 Ccf  +  Da  4-  Ec  =  0 
2  Abe  4  B(bf  +  ed)  4-  2 Cdf  +  Db  +  Ed  =  0 
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These  relations  underdetermine  the  coefficients  a,  b,  c,  d,  e,  /.  Unless  the  transformation  (5.5) 
is  degenerate  we  can  set  b  =  0  and  c  =  0  to  get 

2>4e  +  Bf  =  0 

Be  +  2Cf  +  E  =  0  ^ 

to  determine  e  and  /.  The  coefficients  a  and  d  are  left  undetermined,  and  for  convenience 
they  can  be  set  equal  to  one.  Then  (5.2)  becomes 

A £2  +  Bfr  +  Cy2  +  G  >  0  (5.8) 


where 

G  =  Ae2  +  Bef  +  Cf2  +  De  +  Ef+F  (5.9) 

The  case  where  (5.7)  does  not  provide  valid  values  for  e  and  /  requires  special  treatment 
that  is  not  described  here. 


The  matrix  of  the  quadratic  form  (5.8)  is 


(5.10) 


Let  V  be  the  matrix  whose  columns  are  the  normalized  eigenvectors  of  X.  Then  the  trans¬ 
formation 

fe)=v(0'  (9=vTte)  (5-n) 

gets  rid  of  the  product  term  and  gives 


Hx2  +  Jy2  +  G  >  0 

where 

H  =  AV&  +  BVnVn  +  CV ?2 

j  =  av2\  +  bv21v22  +  cv?2 

Consider  the  toriodal  region.  This  is  the  region  inaccessible  to  the  threat  because  it 
cannot  turn  sharp  enough  to  enter.  This  toroidal  region  is  determined  by  given  values  of  r, 
f,  and  bm.  The  axis  of  the  torus  is  a  circle  of  radius  bm  lying  on  the  plane  perpendicular  to  r 
with  the  point  r  at  its  center.  Each  point  on  this  axial  circle  is  the  center  of  another  circle  of 
radius  bm  that  lies  on  the  surface  of  the  torus  in  the  plane  perpendicular  to  the  axial  circle. 
These  circles  generate  the  torus. 

If  ra  is  a  point  on  the  axial  circle  of  the  torus,  then  the  relations 

I ra  -r\  =  bm 
r  •  (fa  -  r)  =  0 
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must  be  satisfied.  It  is  convenient  to  define  two  special  points  (fao  and  ra\ )  lying  on  the 
axial  circle.  Let  rao  lie  at  the  top  of  the  axial  circle,  and  ra\  lie  on  the  axial  circle  ninety 
degrees  clockwise  from  r„o  relative  to  the  direction  of  r.  The  point  ra o  is  represented  by 


ra  0  =  r  +  ar  +  /3k 

(5.15) 

and  must  satisfy  (5.14).  Accordingly,  the  constants  a  and  /3  must  satisfy 

V2a  +  w(3  =  0 

V2a2  +  2  wajd  +  /?2  =  b2m 

(5.16) 

The  point  ra\  is  represented  by 

r  a\  =  r  +  ai  +  /3j 

(5.17) 

redefining  a  and  (3,  which  are  now  determined  by  the  relations 

ua  +  v(3  =  0 

«2  +  /?2  =  & 

(5.18) 

It  is  convient  to  let  7  be  a  parameter  that  selects  one  of  the  points  ra  on 
and  we  thus  set 

the  axial  circle, 

ra  =  r  +  (foO  —  r)  cos  7  +  (f0i  —  r)sin  7 

(5.19) 

Then  every  7  €  [0,27r)  generates  a  circle  on  the  toroidal  surface.  Each  circle  can  intersect 
the  sea  surface  at  zero  or  two  points.  If  r7  is  one  of  these  points,  there  is  the  relation 


|r7-ra|  =  bm 

rap  •  (fy  -  ra)  =  0  (5.20) 

fr  =  xbi  +  ybj 

(defining  xb  and  yb),  where  fap  is  a  unit  vector  at  ra  along  the  axial  circle.  It  can  be  found 
by  taking  the  derivative  of  (5.19),  giving 

Tap  =  [  {rao  -  r)sin  7  +  (rai  -  r)cos  7 ]/bm  (5.21) 

Figure  5.2  illustrates  the  RA  in  a  typical  case.  The  expected  threat  position  is  indicated 
by  the  star  (between  the  two  circles),  and  its  velocity  is  directed  along  the  dashed  line, 
which  is  the  axis  of  revolution  of  the  figure.  The  HE  restriction  defines  the  conical  surface 
represented  by  the  two  lines  that  would  intersect  at  the  threat  position.  The  toroidal  region, 
represented  by  the  two  circles,  is  inaccessible  to  the  threat  because  of  limitation  on  its  lateral 
acceleration  magnitude  (the  “g”  limitation). 

Figure  5.3  illustrates  the  intersection  of  the  RA  of  Figure  5.2  and  the  sea  surface  for  a 
typical  case  where  the  threat  is  initially  diving.  The  region  of  the  sea  surface  in  the  torus 
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is  shaded.  The  expected  threat  position  is  indicated  by  the  star.  The  projection  of  the 
axial  circle  of  the  torus  onto  the  sea  surface  is  dashed.  The  array  of  dots  lying  just  inside 
the  dashed  curve  are  those  points  on  the  sea  surface  where  the  simplified  ASCM  trajectory 
would  arrive  at  a  point  of  tangency.  The  dots  are  at  equal  increments  of  the  turning  rate  in 
the  hypothetical  simplified  trajectory.  The  outer  solid  curve  indicates  the  boundary  of  the 
conical  region.  The  RA  consists  of  the  unshaded  region  inside  this  solid  curve.  However, 
this  RA  includes  a  subregion  reachable  only  if  the  threat  is  able  to  fly  under  water.  The 
simplified  trajectories  reaching  this  subregion  would  lie  under  the  torus.  The  boundary  of 
this  “wet  trajectory”  subregion  is  not  indicated. 

For  a  given  threat  it  is  interesting  to  consider  those  regions  of  the  sea  surface  where  the 
threat  can  and  cannot  fly.  For  a  given  threat  state  and  the  given  “g”  and  HE  constraints, 
the  preceding  sections  give  rules  that  can  indicate  whether  a  given  point  on  the  sea  surface 
can  be  reached. 

Figure  5.4  shows  several  regions  on  the  plane  of  the  sea  surface  for  a  threat  with  an  initial 
altitude  of  1600  yds,  horizontal  velocity  of  500  yd/sec,  and  vertical  velocity  of  -600  yd/sec. 
The  dive  angle  is  arctan  (600/500)  or  approximately  70  deg.  The  maximum  allowed  lateral 
acceleration  is  taken  to  be  15.8  “g”  units,  giving  a  minimum  allowed  arc  radius  of  3600  yd. 
The  maximum  allowed  HE  (<f>m)  is  taken  as  50  deg.  The  asterisk  shows  the  projection  of  the 
initial  threat  position  onto  the  sea  surface.  The  solid  circle  is  the  “Tangency  Circle”  of  the 
seaskimming  trajectories. 

Regions  3  and  6  are  inaccessible  because  of  the  “g”  constraint,  and  together  represent  the 
intersection  of  a  torus  with  the  sea  surface.  Regions  6  and  7  are  inaccessible  because  of  the 
HE  constraint,  and  together  are  bounded  by  a  conic  section.  Regions  1  and  2  are  accessible 
by  the  “Simple  Trajectories”  described  in  Chapter  3.  Regions  2  and  5  are  accessible  by  the 
“Seaskimming”  trajectories  described  in  Chapter  4.  The  region  accessible  under  the  stated 
assumptions  and  constraints  therefore  consists  of  the  union  of  regions  1,  2,  and  5. 

Regions  4  and  5  would  be  accessible  by  a  “Simple  Trajectory”  if  flight  under  water  were 
allowed.  Region  7  is  accessible  under  the  “g”  constraint,  but  violates  the  HE  constraint  and 
the  condition  that  flight  must  occur  only  in  air  (the  “dry”  constraint). 

Figure  5.5  represents  the  same  conditions  as  Figure  5.4  except  that  <j)m  is  50  deg.  The 
regions  identified  in  Figure  5.4  are  present  in  Figure  5.5,  but  now  there  is  a  new  region. 
Region  8  is  accessible  under  the  “g”  constraint  and  satisfies  the  “dry”  condition  but  violates 
the  HE  constraint. 
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Figure  5.1.  Region  satisfying  the  HE  constraint:  the  cone  radially  symmetric  about  the 
velocity  vector  V  indicated  by  shading. 
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Figure  5.2.  The  region  satisfying  both  the  HE  and  the  “g”  constraints  is  radially  symmetric 
about  the  velocity  vector.  The  shadowing  shows  a  section  of  this  region  about  the  axis  of 
symmetry. 


Figure  5.3.  The  Accessible  Region  on  the  Sea  Surface. 
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Figure  5.5.  Example  of  Accessible  and  Inaccesible  Regions  on  the  Sea  Surface  (with  arbitrary  X  1 0 


NSWCDD/TR-00/46 


CHAPTER  6 

ASCM  DRAG  LOSS 

The  simplified  threat  trajectories  being  used  for  threat  evaluation  consist  of  CV  and 
CTR  segments.  The  drag  and  lift  are  approximated  by 

Fd  =  kV?(Cd  o  +  KC\ ),  FL  =  kV?Cl  (6.1) 

(see,  for  example,  Etkin,  1972)  where  k  is  \pS.  The  balance  offerees  on  the  threat  body  is 
represented  by 

mr  =  FD  +  FL  +  fh  +  mg  (6.2) 

where 

FD  =  -FD?fVtl  FL  =  FLu ,  Th  =  Th?/Vt,  g  =  —gk  (6.3) 

The  lift  is  in  the  direction  of  the  unit  vector,  u ,  which  at  the  moment  is  undetermined. 
However,  it  must  be  perpendicular  to  r  in  addition  to  having  unit  length,  and  therefore  is 
specified  by  only  one  additional  constant. 

The  thrust,  T*,  needed  to  keep  the  threat  speed  constant  is  just  an  artifact  of  the 
simplified  dynamics  and  has  no  meaningful  importance.  It  is  approximated  by  taking  it  in 
the  direction  of  the  threat  velocity  instead  of  along  the  body  axis.  The  purpose  of  this  term 
is  to  avoid  the  necessity  of  considering  the  body  orientation.  The  effect  of  thrust  is  neglected 
because  in  the  first  place  it  is  unknowable,  and  secondly  because  the  gain  in  energy  resulting 
from  any  actual  powered  flight  incurs  a  corresponding  cost  to  the  threat. 

It  may  seem  plausible  to  omit  gravity  in  (6.2)  since  it  would  seem  to  contribute  equally 
to  a  gain  in  kinetic  energy  and  a  loss  in  potential  energy  as  the  threat  descends.  However, 
these  two  effects  do  not  cancel.  Gravity  requires  a  compensating  lift,  which  is  obtained  by  an 
appropriate  angle  of  attack  that  also  affects  the  drag.  Accordingly,  gravity  must  be  included 
in  (6.2).  Furthermore,  it  is  not  necessary  to  account  for  the  gravitational  loss  in  potential 
energy  because  this  loss  is  independent  of  the  trajectory  shape. 

The  problem  is  to  determine  the  total  drag  loss  from  a  threat  trajectory  with  given  r(t) 
and  r(t).  The  drag  as  well  as  the  other  quantities  can  be  found  from  (6.2)  and  (6.3)  by  scalar 
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multiplication  of  (6.2)  successively  by  r,  r,  u,  and  k  to  get 

0  =  -kV?(Cdo  +  KC\ )  +  ThVt  -  mgvz  (6.4a) 

m|r|2  =  kV2Cl(u  •  r)  —  mgaz  (6.46) 

m(u  •  r)  =  kV2Cl  —  mguz  (6.4c) 

maz  =  -KVt(Cm  +  KCl)vz  +  kV2Cluz  +  ThVt~lvz  -  mg  (6.4d) 


where  the  vz,  az,  and  uz  represent  the  z  components  of  r,  r,  and  u.  The  quantities  not  yet 
determined  are  T/,,  uz,  (u-r),  and  Cl.  Only  three  of  the  four  equations  (6.4)  are  independent, 
but  also  only  three  of  the  four  undetermined  quantities  are  independent. 


6.1  The  CV  Segment 

The  acceleration  is  set  equal  to  zero  and  only  three  of  the  equations  (6.4)  are  meaningful. 
Because  of  the  loss  of  one  of  the  equations  (6.4)  the  CV  case  is  not  a  special  case  of  CTR 
and  is  treated  separately.  The  quantities  7\  and  (k  •  u)  are  eliminated  from  (6.4)  to  obtain 

r2  »ng[l  -  {ozJVtf\  - 

L  2nKVtvz  +  K  J 


The  drag  is  then 


Fd  =  K  V?(Cd*  +  KCl) 


with  Cl  as  given  by  (6.5).  The  quantity  (6.6)  is  constant  if  variation  in  air  density  is 
neglected.  The  total  drag  loss  is  just  Fj)  times  the  length  of  the  segment. 


6.2  The  CTR  Segment 

The  quantities  (u  •  r),  7),,  and  (k  •  u)  are  eliminated  from  (6.4)  to  obtain 

c^^?(v'-4+2s6V+sV[i-|])  <6-7) 
A  real  solution  is  obtained  if  (6.7)  is  positive.  Note  that 

sign  {Cl)  =  sign  j  Vt4  |\  +  +  g2b2  \  -  J  (6.8) 

where  a  is  the  acceleration  magnitude.  The  second  bracket  is  always  non-negative.  While 
the  vertical  component  a2  is  often  negative,  there  is  a  wide  range  of  circumstances  where 
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the  first  bracket  is  non-negative,  or  if  it  is  negative  it  could  be  dominated  by  the  positive 
second  bracket.  However,  the  exact  criterion  for  a  real  solution  seems  complicated. 

The  drag  is  then  given  by 

Fd  =  nV?CD o  +  (v*  +  2 gb\,  +  1  -  ^  )  (6.9) 

The  energy  lost  during  the  entire  trajectory  is 


Fj)  ds 


(6.10) 


where  ds  =  d6/b  is  the  incremental  path  length  and  the  integration  is  over  the  entire  CTR 


In  a  CTR  segment  the  threat  position  is  given  by  (3.5),  from  which  the  velocity  and 
acceleration  are  seen  to  be 


r  =  rocos  ut  +  u>(rc  —  ro  )sinu;< 
r  —  —  wrosin  u>t  4-  oj2{rc  —  ro) cos  ujt 


(6.11) 


The  total  loss  is  obtained  by  the  integral  (6.10)  over  0  from  0  to  Of .  Neglecting  the  variation 
in  air  density  (keeping  k  constant),  the  integrals  of  the  time- varying  parts  of  (6.9)  are 


rtt 

I  azdt  =  v2o(cosu )tt  —  1)  -f  uj{zc  —  2o)sin  ujtt 
Jo 


(6.12) 


r*t 

I  v2dt  =  l[v20  +  lo2(zc  —  zq)2]U  +  |u>_1  [v20  -  ui2(zc  —  zo)2]sin  2ojtt 

V  0 

+  jv*o(*c  —  zo)(l  —  cos  2 Ultt) 


(6.13) 


NSWCDD/TR-00/46 


CHAPTER  7 

THE  TBM  TRAJECTORY 


7.1  Introduction 

ASCMs  are  evaluated  by  estimating  an  Objective  Probability  for  each  threat-objective 
pair  by  constructing  a  simplified  trajectory  from  the  initial  threat  position  to  the  objective. 
Many  alternative  simplified  trajectories  were  possible  because  ASCMs  can  maneuver  by 
using  aerodynamic  forces.  On  the  other  hand,  our  treatment  considers  a  TBM  as  incapable 
of  maneuver.  Its  impact  point  is  therefore  predetermined  by  its  initial  threat  state  and  its 
physical  characteristics  such  as  drag  coefficient,  etc.  The  uncertainty  of  its  impact  point 
arises  from  uncertainties  in  these  physical  characteristics  and  errors  in  the  estimation  of  its 
initial  state. 

The  length  of  TBM  trajectories  is  not  negligible  compared  to  the  Earth  radius,  and 
consequently  it  is  advisable  to  use  geographic  rather  than  Cartesian  coordinates.  The  time 
interval  between  establishment  of  firm  track  and  impact  is  long  enough  that  earth  rotation 
needs  to  be  considered.  This  could  be  accomplished  by  bringing  the  Coriolis  force  into  the 
dynamics.  Here,  we  use  the  simpler  approach  of  using  nonrotating  coordinates  and  allowing 
the  earth  to  rotate  under  the  trajectory  during  the  time  in  flight.  To  do  this,  the  initial  E 
velocity  component  is  adjusted  for  Earth  rotation  to  obtain  the  threat  velocity  relative  to 
an  inertial  (nonrotating)  system. 

Assignment  of  Objective  Probabilities  is  effected  by  determining  an  impact  “footprint” 
on  the  Earth’s  surface.  This  footprint  is  the  two-dimensional  PDF  of  the  impact  point.  To 
simplify  the  calculations,  the  expected  impact  point  is  determined  by  the  expected  initial 
threat  state  along  with  the  best  estimates  of  the  threat’s  physical  characteristics.  The 
footprint  is  assumed  to  be  a  Gaussian  function  whose  covariance  is  determined  by  several 
simplifying  assumptions. 

Oblateness  of  the  Earth  is  neglected  and  we  use  the  right-handed  spherical  geographic 
coordinates  (<j>,  L,  r )  representing  East  Longitude,  North  Latitude,  and  radial  distance.  The 
threat  velocity  is  expressed  by  components  (u,  v,  w)  having  the  dimensions  of  L/T  in  the  three 
coordinate  directions.  Note  that  these  components  are  neither  covariant  nor  contravariant. 
Also,  it  is  understood  that  u  has  been  adjusted  for  earth  rotation.  That  is,  (u,  v,  w )  are  the 
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components  in  the  celestial  (inertial)  system.  The  only  difference  between  the  celestial  and 
the  geographic  coordinates  is  in  the  longitudinal  or  E  components,  which  are  related  by 

(j>  =  <f>g  +  ilt,  u  =  ug  + ail  cos  L  (7.1) 

where  the  subscript  “g”  denotes  the  geographic  coordinate  or  component,  and  fl  is  the 
angular  speed  of  the  Earth  rotation,  equal  to  2ir  per  sidereal  day.  The  celestial  longitude 
<f>  is  referred  to  a  datum  fixed  on  the  celestial  sphere.  At  the  initial  time  the  celestial  and 
geographic  coordinate  systems  coincide. 

The  complication  of  atmospheric  reentry  on  a  rotating,  spherical  Earth,  is  avoided  by 
adjusting  a  Keplerian  orbit.  Instead  of  introducing  the  Coriolis  force  into  the  dynamics, 
Earth  rotation  is  considered  by  carrying  out  the  calculations  in  the  celestial  (inertial)  frame. 
The  effect  of  aerodynamic  drag  on  a  planar  trajectory  with  a  flat  Earth  is  used  to  adjust  the 
impact  point  and  flight  time. 

The  following  steps  are  performed: 

•  Use  the  expected  initial  threat  state  relative  to  an  inertial  system.  This  is  done  by 
just  adding  ail  cos  L  to  the  Eastward  velocity  component  relative  to  the  rotating  Earth. 
Calculate  the  arc  length,  0C,  and  flight  time,  T,  to  the  impact  point  on  the  surface  of  the 
spherical  Earth  assuming  a  Keplerian  orbit.  Drag  and  Earth  rotation  are  not  considered. 

•  Carry  out  the  same  calculation  in  Cartesian  Coordinates  on  a  plane  with  a  flat  Earth. 
Find  the  distance  Xq  and  flight  time  To  to  the  impact  point  assuming  a  parabolic  tra¬ 
jectory.  Drag  and  Earth  rotation  are  again  neglected. 

•  Carry  out  the  same  calculation  in  Cartesian  Coordinates  on  a  plane  with  a  flat  earth,  but 
this  time  consider  the  drag  force  in  an  exponential  atmosphere.  Integrate  the  dynamic 
equations  numerically  to  find  the  distance  Xa  and  flight  time  T„  to  the  impact  point. 
The  lift  force  is  considered  to  be  zero. 

•  The  quantities  Xa  —  Xq  and  Ta  —  To  provide  a  correction  for  drag  that  would  be  cum¬ 
bersome  to  obtain  in  spherical  coordinates.  Find  the  corrected  arc  distance  and  time  of 
flight, 

e  =  ec+a-\ X'-x„) 

T  =  TC  +  (T.-T0) 

•  Calculate  the  impact  point  in  the  Celestial  Coordinate  system  based  on  the  expected 
initial  threat  state  relative  to  an  inertial  system,  using  the  corrected  flight  distance  0 
and  flight  time  T. 
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•  Get  the  (rotating)  geographic  coordinates  of  the  impact  point.  This  is  done  by  just 
subtracting  aflT cos  Li  from  the  Celestial  East  longitude  of  the  impact  point. 

•  Use  Cartesian  coordinates  on  a  nonrotating,  flat  Earth  and  parabolic  orbits  to  calculate 
the  (2  by  2)  covariance  matrix  of  the  footprint  PDF  on  the  Earth  surface  from  the  (6 
by  6)  covariance  matrix  of  the  initial  threat-state  errors.  Uncertainties  in  knowledge  of 
the  threat  characteristics  (drag  coefficient,  etc.)  are  included  to  increase  the  footprint 
covariance. 

•  Assume  a  Gaussian  footprint  on  the  Earth's  surface  having  the  expected  position  and 
covariance  as  obtained  above. 

The  last  two  steps  are  described  in  Chapter  8. 


7.2  The  Keplerian  Orbit 

Refer  to  Figure  7.1.  With  the  simplifications  described  above  the  arc  length  6  between 
the  initial  threat  geographic  position  ( <f>0,L0 )  and  the  impact  point  (&,  Li)  depends  only  on 
the  initial  threat  speed  Vo  and  flight-path  angle  70.  This  relation  is  given  by  Regan  and 
Anandakrishnan  (1993,  equation  (6.34))  as 

^  _  1  —  cos  0 

(ro/ Re)cos27o  —  cos(#  +  7o)cos  70 

where  ro  =  a  -f  zo  is  the  initial  threat  radial  distance,  a  is  the  Earth  radius,  Re  =  a  +  2;  is 
the  radial  distance  of  the  impact  point,  and 


A  = 


Yti 


ft  =  goa2 


(7.4) 


where  go  is  the  acceleration  of  gravity  at  the  Earth  surface.  The  parameter  A  determines  the 
trajectory  type.  For  the  suborbital  trajectories  of  interest,  A  €  (0, 1).  Solving  this  equation 
for  0  by  representing  the  cosine  in  terms  of  the  sine  leads  to  the  equation 


A  sin2#  +  2?  sin#  +  C  =  0  (7-5) 

where 

A  —  A2cos27o  —  2A  cos27q  +  1 

B  =  -2A  cos  70  sin  70(1  —  A<rcos27o)  (7.6) 

C  =  A  (a  —  1)cos27q  [(er  +  1)A  eos27o  —  2] 

where  a  is  ro  /  Re .  It  should  be  noted  that  there  are  cases  where  the  trajectory  may  intersect 
the  Earth  surface  behind  the  initial  threat  position,  or  may  not  intersect  the  Earth  surface 
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at  all.  These  cases  are  rejected,  and  the  solution  of  interest  is  the  smallest  positive  value  0 
obtained. 


The  unadjusted  time  of  flight  is  given  by  Regan  and  Anandakrishnan  (1993,  equation 
(6.41b)),  and  is 


Tf  = 


r  o 


{ 


tan  70(1  —  cos  0)  +  (1  —  A)sin  0 


V0COS70  i  (2  -  A){[(1  -  cos0)/(Acos27o)]  +  [cos(7o  +  0)/cos7o]} 


} 


+ 


2ro 


V0A[(2/A)  -  I]3/2 


arctan 


[(2/A)  -  l]1/2 


cos  70  ctn(0/2)  —  sin  70 


(7.7) 


7.3  Correction  for  Drag 


It  is  assumed  that  the  reentry  vehicle  is  subjected  to  drag  but  not  to  any  lift  force.  The 
simplified  aerodynamic  equations  are 

dx  t/ 

—  =  V  cos  7 

dt 


dz 

—  =  V  sin  7 

dt  ' 

^  =  -2V[8(e~z/h  +  g  sin  7] 

T/-1 

-  =  -gV  cos  7 


(7.8) 


where 


<  =  v 2,  6  = 


pgCpS 

2m 


(7.9) 


with  po  being  the  air  density  at  the  Earth’s  surface,  Cj)  the  drag  coefficient,  S  a  cross- 
sectional  area,  and  m  the  threat’s  mass.  These  equations,  which  have  been  adapted  from 
Regan  and  Anandakrishnan  (1993),  are  integrated  numerically  to  obtain  the  time  of  flight,  Ta, 
and  horizontal  distance  traversed,  Xa,  by  a  body  falling  through  the  exponential  atmosphere. 


Modification  of  the  Keplerian  trajectory  is  made  by  comparing  the  cases  of  falling  through 
air,  as  in  the  previous  paragraph,  and  falling  through  a  vacuum.  For  the  latter  case  the  time 
of  flight  and  horizontal  distance  traversed  are  given  by 


To 


zo  +  yjz\  +  2p(z0  -  Zi) 


Xq  =  xqTq 


(7.10) 


to  provide  the  quantities  needed  in  (7.2)  for  this  modification.  An  example  is  shown  in 
Figure  7.2. 
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7.4  Impact  Point 

The  foregoing  formulas  are  for  determining  the  impact  point  for  a  threat  whose  initial 
position  and  velocity  are  given.  Consider  the  spherical  triangle  NOS  where  N  is  the  North 
Pole,  O  is  the  subthreat  point,  and  S  is  the  subimpact  point,  where  wsub”  means  projected 
radially  onto  a  sphere  of  radius  a  (see  Figure  7.3).  From  the  previous  definitions  it  is  seen 
that  2"  —  L0  is  the  arc  NO,  0  is  the  arc  OS,  and  ^  —  L,  is  the  arc  NS.  Let  the  angle  NOS  be 
denoted  by  a,  where 

tan  at  =  -uo/v0  (7,11) 

Then  the  celestial  coordinates  (fa,  Li)  of  the  impact  points  are  given  by 

sin  Li  =  sin  Lq  cos  8  +  cos  Lq  sin  0  cos  a  (7-12) 

and 

.  ,  ,  ,  .  sin  8  sin  a 
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KEPLERIAN  ORBIT 


Figure  7.1.  The  Keplerian  orbit.  The  flight  plane  passes  through  the  threat  initial 
position,  the  Earth  center,  and  the  initial  velocity  vector.  This  figure  is  viewed  in  an  inertial 
reference  frame,  in  which  the  Earth  is  rotating. 
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•  Figure  7.2.  Example  of  fall  through  exponential  atmosphere,  and  through  vacuum  The 
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NP 


Figure  7.3  Location  of  the  impact  point  “I”  on  the  spherical  Earth,  where  “0”  is  the 
projection  of  the  initial  threat  position  onto  the  Earth  surface  and  “NP”  is  the  North  Pole. 
This  picture  is  as  viewed  in  an  inertial  reference  frame,  with  the  Earth  rotating. 
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Example  TBM  Trajectory  and  Footprint 


East  Longitude  (deg) 


Figure  8.2.  Example  TBM  Trajectory  Projected  onto  Sea  Surface  and  Footprint. 
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CHAPTER  8 

THE  TBM  FOOTPRINT 


8.1  Covariance  of  Impact  Position  Error 

The  geometry  of  Keplerian  trajectories  in  spherical  coordinates  is  more  complicated  than 
necessary  for  the  estimation  of  the  impact  position  error.  For  this  estimate  we  consider  a 
parabolic  trajectory  on  a  plane  Earth  in  a  right-handed  Cartesian  system  (x,y,z)  in  the 
direction  of  the  local  (<f>,L,r)  (with  z  directed  upward),  and  corresponding  velocity  compo¬ 
nents  ( x,y ,  z).  Also,  the  acceleration  of  gravity,  g,  is  assumed  constant.  The  flight  time  and 
the  impact  point  at  elevation  Z{  are  given  by 


Tf  =  g 


-  „-i 


io  +  yj. z\  +  2g(z0  -  zi) 


x%  =  x0  +  xqTf,  yi  =  yo  +  yoTp  (8.1) 


Even  if  the  initial  state  error  is  Gaussian  the  joint  PDF  of  the  position  error  of  the  impact 
point  is  a  complicated  function  because  of  the  nonlinearity  of  (8.1).  To  simplify  the  analysis, 
the  initial  state  errors  are  now  taken  to  be  infinitessimal  quantities.  To  relate  the  initial 
state  errors  with  the  impact  position  errors  we  need  the  (2  by  6)  Jacobian  matrix  J,  where 


=  (  fl(g»’»y«)  \ 

\d(xo,  yo,  zo,  io,  yo,  h) ) 


( 


1 

0 


io  [io  +  2y(2o  -  *»)]  1/2 

r  - i - 

io  +  yig  +  2y(20-2,) 


0 

1 


\ 


yo  [io  +  2y(zo  -  zi)] 
0 


0  9  1  Jio  +  yjz%+  2g(z0  -  Zi) 

^y-^o  (l  +  i0[io  +  2g(zo  -  2,)]“1/2)  g~lm  (l  +  io[io  +  2 g(z0  -  2,)]-1/2) 


(8.2) 


} 


The  (2  by  2)  covariance  of  the  impact  position  error  Cs  due  to  uncertainty  in  the  threat 
state  is  related  to  the  (6  by  6)  covariance  of  the  initial  threat  state  error  Co  by 

Cs  =  JC0JT 
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8.2  Uncertainty  in  the  Drag  Coefficient 

If  the  threat  classification  is  unknown,  there  may  be  appreciable  uncertainty  in  the 
appropriate  value  of  the  drag  parameter  8.  This  uncertainty  will  widen  the  extent  of  the  TBM 
footprint  and  needs  to  be  taken  into  account.  There  are  ways  to  estimate  the  drag  and  its 
uncertainty  from  previous  state  estimates.  A  few  such  methods  are  outlined  in  Appendix  A. 
The  following  is  a  simplified  procedure  for  estimating  the  TBM  footprint  variance  assuming 
that  the  expected  value  and  standard  error  of  5  are  given. 

Let  8  and  ag  be  the  mean  and  standard  deviation  of  the  hypothetical  PDF  of  8  represent¬ 
ing  our  ignorance  of  its  correct  value.  Carry  out  the  numerical  integration  of  the  equations 
(7.8)  for  two  values  of  8  and  denote  the  results  as  follows: 


Tp\ i  for  8  =  8  —  og 

Tp,  0  for  8  =  8 


(8.4) 


Denote  the  geographic  coordinates  of  the  two  impact  points  by  (0i,  L\)  and  (0,  L).  The  un¬ 
certainty  in  local  Cartesian  coordinates  due  to  the  error  in  8  is  represented  by  the  covariance 
matrix 

C„  =  a2(,  <*r,?)26in'i  ,  (■/>!- fm  -  L)sin  L\  (8.5) 

V  (01  —  0)(^i  ~  L) sin  L  (L\  —  L)  J 

The  total  error  covariance  in  the  impact  point  is  given  by 


Ci  =  CS  +  Cd  (8.6) 

The  impact  “footprint”  resulting  from  the  initial  state  vector  and  error  covariance  is  taken  to 
be  a  Gaussian  PDF  with  mean  given  by  the  geographic  impact  location  ((7.12)  and  (7.13)) 
of  the  expected  threat  state,  and  covariance  given  by  (8.6).  Some  examples  are  shown  in 
Figures  8.1  and  8.2.  When  the  TBM  footprint  overlaps  a  defended  area,  as  it  always  does  to 
some  extent,  as  in  Figure  8.3,  the  probability  that  the  TBM  will  impact  within  the  defended 
area  is  the  integral  of  the  footprint  PDF  over  the  defended  area,  and  represented  by  Pkm,  as 
discussed  in  Section  9.3. 
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CHAPTER  9 

DETERMINATION  OF  EVALUATION  PARAMETERS 

9.1  Objective  Probabilities 

The  immediate  problem  is  to  determine  the  objective  of  each  threat  in  a  battle  situation. 
An  objective  is  any  surface  unit  or  other  defended  point /region  that  may  be  the  objective 
of  an  enemy  attack.  In  situations  where  there  are  a  large  number  of  defended  positions 
that  might  be  under  attack,  as  in  the  case  of  a  task  force,  there  is  uncertainty  as  to  which 
position  is  the  objective  of  any  given  threat.  This  suggests  the  assignment  of  an  “Objective 
Probability”  to  each  threat-objective  pairing.  Let  Pkj  denote  the  probability  that  objective 
3  is  the  target  of  threat  k.  Then  for  a  threat  capable  of  reaching  any  of  the  objectives 

Y  Pkj  =  1  for  all  k  (9.1) 

j 

Otherwise 

Pkj  =  0  for  all  j  (9.2) 

for  those  threats  k  incapable  of  reaching  any  objective.  Then  we  can  use  a  Markov  update 
model  to  assess  probability  that  a  state  will  be  reached  by  a  threat  k. 

9.2  ASCMs 

The  targets  of  ASCMs  are  assumed  to  be  defended  points,  each  having  a  specified  po¬ 
sition.  Some  of  these  positions  may  change  during  an  engagement.  The  probability  that  a 
given  threat  is  targeting  a  given  objective  would  seem  to  depend  on  several  circumstances: 
(1)  the  perceived  value  of  the  objective,  (2)  the  flight  time  to  the  objective,  and  (3)  the 
energy  requirement  to  reach  the  objective. 

Each  defended  point  is  assigned  a  Perceived  Value  vj^  and  an  Actual  Value  vja\  The 
value  vjp)  determines  the  probability  of  its  being  targeted,  while  vja)  determines  the  effort 
to  be  expended  in  its  defense.  These  values  are  not  necessarily  equal,  because  the  defender 
may  have  some  knowledge  about  the  objective  that  the  attacker  lacks. 

The  flight  time  and  energy  requirement  depend  on  the  unknown  threat  trajectory.  Ways 
to  hypothesize  the  trajectory  of  an  ASCM  and  arrive  at  plausible  values  of  the  flight  time 
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and  energy  requirement  are  discussed  in  Chapters  3,4,  and  6.  Let  Tkj  be  the  ETA  flight  time 
of  threat  k  to  objective,;'.  Similarly,  let  L^j  be  the  drag  loss  incurred  by  threat  k  in  its  ETA 
flight  to  objective  j.  These  values  are  used  to  calculate  a  preliminary  objective  probability 
Pkj.  It  seems  reasonable  to  assume  that  the  probability  of  targeting  an  objective  increases 
both  with  flight  time  and  drag  loss.  Accordingly,  the  formula 


9.3  TBMs 

For  each  threat  its  initial  state  estimate  determines  a  “TBM  footprint”  in  terms  of  the 
PDF  of  its  impact  point,  as  described  in  Chapters  7  and  8.  The  assignment  of  values  of  P^j  to 
the  defended  positions  is  different  for  defended  points  and  defended  areas.  Consider  first  the 
defended  points.  Each  defended  point  lying  within  this  footprint  is  assigned  a  preliminary 
value  of  Phj  in  proportion  to  the  footprint  density  at  its  position. 

Each  defended  region  is  assigned  a  Perceived  Value  per  Unit  Area  wjp^  and  an  Actual 
Value  per  Unit  Area  wja\  The  reason  for  defining  the  values  of  defended  areas  differently 
from  those  for  defended  points  is  explained  below.  Defended  areas  have  the  complication 
that  the  radius  of  destruction  about  the  impact  point  depends  on  the  yield  of  the  threat, 
which  may  not  be  known.  The  simplified  procedure  described  here  considers  this  effect  only 
through  values  of  the  “Threat  Menace”  Af*,  which  are  assumed  specified  along  with  other 
features  of  threat  classification  in  advance  of  the  “Threat  Evaluation”  being  described  here. 

Each  defended  region  lying  within  or  overlapping  the  footprint  is  assigned  a  preliminary 
value  of  Pkj  in  proportion  to  the  probability  that  the  TBM  will  impact  within  the  area.  This 
probability  is  a  surface  integral  of  the  footprint  PDF  over  the  defended  area. 

The  values  obtained  need  to  be  normalized,  noting  that  the  dimensions  of  the  preliminary 
values  for  the  defended  points  are  different  than  those  for  the  defended  areas.  However,  the 
product  of  the  preliminary  values  P^j  with  the  corresponding  wjp^  have  the  same  dimensions 
as  the  preliminary  P^j  values  for  the  defended  points.  Accordingly,  values  P \j  =  wjp^  Pf.j 
are  defined. 

9.4  Normalization 

Now  that  the  dimensions  of  the  are  the  same  for  defended  points  and  defended  areas, 
these  are  all  normalized  together  so  that  their  sum  over  all  defended  positions  is  unity  for 
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each  threat.  That  is, 


Pt,  =  hi 


for  all  k 


(9.4) 


9*5  Time-Ordered  Evaluation  Parameters 

It  would  be  gieat  if  the  threats  could  be  adequately  characterized  by  a  single  parameter 
each.  However,  it  is  easy  to  visualize  situations  in  which  such  a  simple  characterization 
fails.  Consider  two  threats,  “A*  and  “B*,  where  threat  “A*  reaches  a  low- valued  asset 
quickly,  while  B  leaches  a  more  valuable  asset  later.  Even  in  this  simplistic  scenario  it  is 
evident  that  the  Weapons  Control  System  (WCS)  needs  information  that  includes  times  of 
arrival.  In  addition,  the  optimal  use  of  defensive  assets  (missiles)  requires  additional  types  of 
information,  such  as  estimates  of  kill  probabilities  for  the  various  threats,  shooting  platforms, 
launch  times,  etc,,  not  to  mention  the  scheduling  of  illuminators  and  magazine  limitations. 

The  problem  addressed  here  is  to  provide  an  evaluation  of  the  threats  based  solely  on 
the  current  estimates  of  their  states  and  properties,  along  with  the  positions  and  values 
of  the  defended  assets.  The  evaluation  parameters  provided  by  this  study  are  insufficient 
for  prioritization  of  the  threats,  which  is  expected  to  be  carried  out  subsequently.  Defense 
consists  of  the  following  steps: 


THREAT  EVALUATION  PRIORITIZATION  ENGAGEMENT 

Consider  now  the  product  that  “THREAT  EVALUATION*  should  provide  to  “PRIOR¬ 
ITIZATION  .  Two  alternative  characterizations  are  proposed:  (1)  The  “Maximum  Possible 
Value  Loss*  (MPVL),  and  (2)  the  “Maximum  Expected  Value  Loss*  (MEVL). 

The  MPVL,  for  each  threat,  is  a  monotonic  function  of  time  beginning  with  the  current 
time.  At  any  given  time  the  MPVL  is  equal  to  the  value  (F<«))  of  the  most  valuable  asset 
that  the  threat  could  reach  by  means  of  an  ETA  trajectory.  At  this  time  the  threat  is  capable 
of  reaching  a  definite  set  of  assets.  This  set  includes  the  assets  belonging  to  previous  sets. 
Since  the  MPVL  is  the  maximum  F^  in  the  set,  it  monotonically  increases  with  time. 

The  earliest  arrival  time  of  threat  k  to  defended  asset  j  is  Tkj.  For  each  threat  k  these 
values  are  reordered  to  obtain  a  set  of  ordered  arrival  times  7\.j,  where 

Tk i  <  Tk2  <  Tkz  <  ...  (9.5) 

where  every  i  corresponds  to  one  of  the  asset  indices  j  for  the  given  threat  k.  Having  ordered 
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the  arrival  times  to  obtain  the  pairings 


(k,i)  ~  (k,j)  (9.6) 

the  ordered  Actual  Asset  Values  are  obtained  from  the  original  vja^  by  the  pairing 
(9.6).  The  MPVL  values  are  then  obtained  by  taking 

Eki  =  m&x(Ekl,  Em,  ...,  Eki)  (9.7) 

The  MEVL  is  obtained  in  the  same  way,  by  replacing  the  vja^  values  by  vja^Pkj.  Note 
that  the  pairings  (9.6)  and  the  Tki  are  the  same  as  before.  However,  the  Ek{  values  may  be 
different  than  in  the  MPVL  case. 

The  implementation  of  threat  evaluation  is  capable  of  providing  both  the  MPVL  and  the 
MEVL.  Neither  is  voluminous  or  not  very  computer  intensive.  It  is  suggested  that  simulation 
with  various  scenarios  be  performed  in  order  to  decide  which  of  the  two  is  preferable. 

9.6  Algorithm 

The  determination  of  the  Objective  Probabilities  proceeds  in  the  following  order: 

•  Determine  the  positions,  and  the  vjv'  and  vja^  values  of  all  defended  points,  and  the 
locations,  extents,  and  the  and  wja^  values  for  all  defended  areas. 

•  Obtain  from  the  Track  Server  the  estimated  initial  states  and  covariance  matrices  for  all 
active  threats. 

•  Obtain  or  hypothesize  classification  data  on  each  threat.  If  an  ASCM,  specify  its  Mk, 
<j>(™\  and  4m).  If  a  TBM,  specify  its  Mk,8k,  and  A 5k. 

•  For  each  ASCM,  k,  determine  the  simplified  trajectory  to  each  defended  point  j  and  the 
values  Tkj,  Lkj ,  and  Pkj. 

•  For  each  TBM,  k,  determine  its  footprint,  and  the  values  Tkj  and  Pkj  for  each  defended 
point  and  area. 

•  Determine  the  Tk  and  Eki  values. 

9.7  Example 

A  task  force  consists  of  six  units  deployed  on  a  two  by  three  rectangular  array  as  in  Figure 
9.1  (on  an  arbitrary  scale),  with  separation  of  three  kilometers  between  adjacent  units.  A 
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threat  is  initially  observed  at  an  altitude  of  five  kilometers,  at  a  distance  of  five  kilometers 
downrange.  The  threat  is  initially  in  horizontal  flight  and  arrives  in  a  circular  arc,  targeting 
unit  number  five.  The  relative  values  of  the  units  values)  are  (1,5, 2, 2, 4, 2),  respectively. 

The  variation  with  time  of  the  “Time  to  Go”  for  the  ETA  trajectory  to  each  unit  (Tj 
values)  is  displayed  in  Figure  9.2  on  an  arbitrary  scale.  These  values  for  all  units  decrease 
as  the  threat  approaches,  but  finally  only  the  Tj  for  unit  five  (the  targeted  unit)  decreases 
to  zero,  as  the  threat  arrives  after  13  seconds.  The  variation  with  time  of  the  “Objective 
Drag  Loss”  for  the  ETA  trajectory  to  each  unit  (Lj  values)  is  displayed  in  Figure  9.3  on  an 
arbitrary  scale.  These  values  exhibit  about  the  same  trends  as  the  Time  to  Go  values. 

The  variation  with  time  of  the  Objective  Probabilities  for  each  unit  (Pj  values)  is  shown 
in  Figure  9.4  on  an  arbitrary  scale.  These  values  were  computed  according  to  equations  (9.3) 
and  (9.4).  The  defenders  at  first  would  assume  that  it  is  more  probable  that  the  threat  is 
targeting  unit  2  because  of  its  higher  value.  But  after  about  three  seconds,  the  Objective 
Probability  of  unit  5  (the  true  target)  becomes  higher,  as  the  Time  to  Go  and  Objective 
Drag  Loss  to  all  other  units  become  great  enough  to  depress  their  Objective  Probabilities  in 
comparison  to  that  for  unit  5,  the  true  targeted  unit. 

The  evaluation  parameters  MPVL  and  MEVL,  as  described  in  Section  9.5,  are  illustrated 
in  another  example  described  in  Appendix  B. 
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Example  I 


Figure  9.1.  Scenario  with  six  defended  points  and  one  ASCM  threat.  Defended  Point 
#5  is  the  threat’s  target,  initially  unbeknownst  to  the  defenders. 
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Time 


Figure  9.3.  Evolution  of  the  ETA  drag  loss  for  the  scenario  of  Figure  9.1. 
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CHAPTER  10 

ORGANIZATION  OF  THE  COMPUTATIONS 

10.1  Background 

The  algorithm  described  in  this  chapter  considers  two  types  of  threat  -  ASCMs  and 
TBMs  -  as  well  as  two  kinds  of  defended  positions  -  Points  (e.g.  ships)  and  Shore  Areas 
(e-9-  cities).  Treatment  of  ASCMs  differs  from  treatment  of  TBMs  in  the  following  respects. 
Each  ASCM  is  considered  potentially  able  to  target  any  of  the  Defended  Points.  These  are 
given  “Objective  Probabilities”  to  indicate  which  are  more  likely  to  be  targeted. 

For  simplicity,  the  Objective  Probabilities  are  based  on  the  expected  state  of  the  ASCM, 
regardless  of  the  errors  in  the  state  estimation.  A  maneuver  is  required  for  the  threat  to  reach 
any  of  the  Defended  Points.  A  hypothesized  simplified  trajectory  is  proposed  for  each  threat- 
objective  pair,  and  used  to  determine  the  time  of  arrival  and  energy  usage  requirement.  The 
Objective  Probability  is  high  for  points  that  can  be  reached  quickly  (early  time  of  arrival) 
and  easily  (low  drag  loss).  The  Objective  Probability  is  zero  for  those  points  impossible  for 
the  threat  to  reach  under  reasonable  assumptions  regarding  its  flight  capabilities. 

In  contrast,  the  trajectory  of  a  TBM  is  considered  determined  by  the  threat  character¬ 
istics  (e.g.,  mass,  drag  coefficient,  etc.)  and  initial  state.  Maneuvers  are  not  considered.  If 
its  target  were  considered  determined  only  by  its  expected  state,  most  likely  it  would  not 
reach  any  of  the  Defended  Points.  These  are  possible  targets  only  because  the  initial  threat 
state  and  its  physical  properties  are  not  precisely  known.  Accordingly,  both  the  covariance 
matrix  of  the  threat-state  estimation  errors  and  uncertainty  in  the  value  of  the  threat’s 
drag  coefficient  are  considered.  These  determine  a  “footprint”,  which  is  the  PDF  of  the 
threat’s  impact  point.  Defended  Points  lying  within  this  footprint  are  assigned  Objective 
Probabilities  related  to  the  footprint  PDF. 

Extended  shore  areas  may  be  targets  of  TBMs,  but  not  of  ASCMs  unless  there  are 
Defended  Points  within  the  area.  The  Defended  Points  are  represented  by  their  positions,  . 
which  may  move  during  the  engagement,  and  by  an  “Actual  Value”  (Va),  and  a  “Perceived 
Value”  (Vp).  The  Perceived  Value  is  used  to  calculate  the  Objective  Probability,  while  the 
Actual  Value  determines  how  much  effort  is  to  be  expended  in  its  defense. 
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A  Defended  Area  can  be  represented  by  a  circle,  ellipse,  or  polygon.  A  TBM  generates 
an  Objective  Probability  determined  by  the  intersection  of  its  footprint  onto  the  Defended 
Area.  There  is  no  structure.  A  real  defended  area  such  as  a  city  would  have  variation  in 
value  from  one  part  to  another.  This  feature  can  be  simulated  by  placing  defended  points 
within  the  area. 

The  Objective  Probabilities  form  a  two-dimensional  matrix  whose  dimensions  equal  the 
number  of  threats  (active  at  a  given  time)  in  one  dimension,  and  the  number  of  Defended 
Points  and  Defended  Areas  in  the  other.  These  values,  along  with  the  Va  values,  constitute 
an  important  input  into  the  evaluation  of  the  threats. 

Engageability  of  the  threat  is  not  considered.  The  number  and  positions  of  the  shooting 
platforms  is  not  included  among  the  data  used.  No  PIPs  are  generated.  The  Objective 
Probabilities  are  estimated  without  regard  to  whether  or  not  the  threat  can  be  engaged,  or 
the  most  favorable  launch  times,  etc. 

10.2  Simplified  Trajectories 

For  each  (ASCM  -  Defended  Point)  pair  a  trajectory  producing  the  soonest  arrival  time 
is  determined.  The  initial  threat-state  error  covariance  is  ignored,  and  the  trajectory  is  based 
entirely  on  the  initial  expected,  state.  The  initial  speed  is  maintained  constant  throughout 
the  flight.  The  trajectory  producing  the  earliest  arrival  time  is  determined.  This  lies  on  the 
plane  determined  by  the  initial  threat  velocity  and  the  Defended  Point.  Initially  there  is  a 
maximum-g  turn  followed  by  a  CV  segment  to  the  Defended  Point.  Cases  where  this  pattern 
would  require  the  threat  to  fly  below  the  sea  surface  are  adjusted  if  possible. 

The  initial  state  estimate  of  a  TBM  including  the  covariance  matrix  generates  a  suite  of 
hypothetical  trajectories,  one  for  each  possible  initial  state.  The  flight  follows  a  Keplerian 
orbit  modified  for  liftless  drag  in  an  exponential  atmosphere. 

10.3  Determination  of  the  Objective  Probabilities 

At  any  given  time,  all  the  threats  that  are  active  at  that  time  are  considered  one  by  one. 
The  ASCMs  are  treated  separately  from  the  TBMs. 

For  an  ASCM,  all  the  Defended  Points  are  considered.  For  each  Defended  Point  the 
drag  loss  along  the  simplified  trajectory  is  approximated  on  the  basis  of  whatever  threat 
information  is  at  hand.  All  incomplete  threat  information  is  hypothesized  in  a  way  that 
minimizes  flight  time  and  drag  loss.  The  Objective  Probability  for  each  Defended  Point  is 
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based  on  these  two  values.  After  all  Objective  Probabilities  have  been  approximated  in  this 
way,  they  are  normalized  so  that  their  sum  over  all  Defended  Points  is  unity  for  each  threat. 

For  a  TBM,  the  first  step  is  to  approximate  its  footprint  on  the  Earth’s  surface.  The 
expected  initial  threat  state  is  used  to  determine  a  single  Keplerian  orbit  to  the  Earth 
surface.  Next,  the  atmospheric  effect  is  approximated  by  using  a  simplified  dynamic  model 
of  a  mass  point  falling  in  a  planar  orbit,  without  lift,  from  the  expected  initial  state  to  the 
Earth  s  surface.  The  resulting  flight  time  and  distance  provide  a  correction  to  the  previously- 
determined  Keplerian  orbit  to  give  a  more  accurate  impact  position  and  flight  time. 

The  initial  state  covariance  is  usually  assumed  to  be  Gaussian.  But  even  with  this 
idealization,  the  complicated  nonlinear  flight  path  makes  it  difficult  to  get  a  closed-form 
footprint  PDF.  Instead,  the  impact  position  covariance  is  determined  by  assuming  small 
errors  and  a  parabolic  planar  trajectory  through  an  airless  medium.  In  addition  to  the 
error  in  impact  position  due  to  error  in  the  estimation  of  the  initial  threat  state,  the  error 
due  to  ignorance  of  the  threat’s  drag  coefficient  is  estimated.  This  is  done  by  making  drag 
corrections  corresponding  to  two  values  of  the  drag  coefficient  differing  by  one  standard 
deviation  of  the  uncertainty.  The  separation  of  the  two  impact  points  provides  an  increment 
to  the  footprint  covariance.  This  impact  covariance  is  centered  on  the  expected  impact 
point  that  was  previously  determined  by  a  somewhat  more  complete  analysis,  and  is  used  to 
generate  a  two-dimensional  Gaussian  PDF  to  represent  the  impact  footprint. 

Defended  Points  lying  within  the  impact  footprint  are  assigned  Objective  Probabilities 
related  to  the  Perceived  Value  and  the  PDF  density  at  that  position.  Defended  Areas 
intersecting  the  footprint  are  assigned  Objective  Probabilities  related  to  their  Perceived 
Values  and  the  integrated  footprint  probability  of  impact  in  the  intersection. 

10.4  Pseudocode  for  Global  Threat  Evaluation 

Initialize  with  n  =  0  (no  threats),  and  m  =  0  (no  defended  positions).  Actions  to  be 
performed  are  one  of  the  following: 

•  Add  a  threat  track. 

•  Drop  a  threat  track. 

•  Update  a  threat  track. 

•  Add  a  defended  point. 

•  Delete  a  defended  point. 

•  Move  a  defended  point. 
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•  Add  a  defended  area. 

•  Delete  a  defended  area. 

•  Calculate  the  Threat  Evaluation  Parameters. 

perform  a  new  action 

update  the  battle  time  to  the  time  of  the  action 
if  new  track  added 
n  —  n  +  1 
state(n)  inputted 
M (n)  inputted 
elseif  track  dropped 
if  kd  <  n 
Mkd  =  Mn 
stat  e(ka)  —  state(n) 

end 

n  =  n  —  1 

elseif  threat  track  ky,  is  updated 
state(fcu)  updated 
Mku  updated 

elseif  a  new  defended  point  added 
m  =  m  +  1 

position,  Vm  \  and  v4°*  inputted 
elseif  a  defended  point  jd  is  deleted 
if  jd  <  m 

position(jd)  =  position(m) 
parameters^)  =  parameters(m) 

end 

m  =  m  —  1 

elseif  defended  point  jm  is  moved 
position(jim)  is  changed 
elseif  a  defended  area  is  added 
m  =  m  +  1 

position,  Wm\  and  inputted 
elseif  defended  area  jd  is  deleted 
if  3d  <  m 
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position^)  =  position(m) 
parameters  (jd)  =  parameters(m) 

end 

to  =  to  —  1 

elseif  Threat  Evaluation  Parameters  are  to  be  calculated 
for  k  =  1  :  n 

if  threat  (fc)=ASCM 
for  j  =  1  :  to 

if  position  (j)  is  a  defended  point 

determine  an  ETA  trajectory  from  k  to  j 
if  the  ETA  trajectory  is  valid 
calculate  Tjy  and  Lkj 

end 

end 

end  (of  j  loop) 
calculate  Pkj 
if  threat  (fc)=TBM 

calculate  the  Keplerian  orbit,  range  and  time 
correct  the  Keplerian  parameters  for  drag 
calculate  the  TBM  footprint 
for  j  =  1  :  to 

if  position  (j)  is  a  defended  point 
calculate  P^j 

elseif  position  (j)  is  a  defended  area 
calculate  P^j  and  hj 

end 

end  (of  j  loop) 
end  (of  if  ASCM  or  TBM) 
calculate  Pfcj.  Tjfc,  and  Eij 
end  (of  k  loop) 
end  (of  action) 

The  organization  of  these  computations  is  shown  schematically  in  Figure  10.1. 
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10.5  Computation  Modules 

The  following  numbers  are  defined: 

na  Number  of  Cruise  Missiles  (ASCMs) 
nt  Number  of  TBMs 
mp  Number  of  Defended  Points 
ma  Number  of  Defended  Areas 
nd  Number  of  Data  Points  in  Trajectories 

derivs  -  a  function  used  during  Runge-Kutta  numerical  integration  (see  Press  et  al,  1986) 
Inputs 

t  time  (s) 

y  (4  by  1)  Vector  of  Dependent  Variables 

delta  Modified  Drag  Coefficient  (1/m) 
hatm  1/e  altitude  of  air  density  (m) 

Outputs 

dydt  (4  by  1)  Vector  of  Derivatives 


dragloss  -  a  function  to  calculate  the  drag  loss 


Inputs 

Cdo 

K 

rOvec 

rsvec 

rOdotvec 

bm 

kappa 

m 

flag 

etarcvec 

etartvec 

etathetat 

Outputs 

etaloss 


Drag  Coefficient  at  Zero  Lift 
A  Coefficient  in  Drag  Relation 
(3  by  1)  Initial  Threat  Position  (m) 

(3  by  1)  Targeted  Point  Position  (m) 

(3  by  1)  Initial  Threat  Velocity  (m/s) 
Minimum  Allowed  Turn  Radius  (m) 

Physical  Parameter  for  Air  and  Vehicle  (kg/m) 
Vehicle  Mass  (m) 

Flag  Indicating  Validity  of  Inputs 

(3  by  1)  Position  of  ETA  Arc  Center  (m) 

(3  by  1)  Position  of  Transition  Point  (m) 

Arc  Length  (rad) 


Drag  Loss  (kg  m2/s2) 
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falltime  -  a  function 

Inputs 

hatm 

delta 

VOcel 

gamqOcel 

zO 

zi 

dt 

epsi 

hi 

hmin 

Outputs 

Tfall 

Xfall 

footprin  -  a  function 

Inputs 

elongO 

nlatO 

rOdotvec 

CO 

zO 

zi 

hatm 

delta 

dt 

epsi 

hi 

hmin 

Outputs 

elongi 

nlati 


called  for  each  TBM  at  desired  time 


Altitude  where  air  density  is  1  /e  of  surface  value  (m) 
Drag  Coefficient  (1/m) 

Initial  Threat  Speed  (m/s) 

Initial  Flight  Angle  (rad) 

Initial  Threat  Altitude  (m) 

Altitude  of  Impact  Area  (m) 

Time  Increment  for  Numerical  Integration  (s) 
Parameter  Needed  in  Runge  Kutta  Integration 
Parameter  Needed  in  Runge  Kutta  Integration 
Parameter  Needed  in  runge  Kutta  Integration 


Time  of  Flight  thru  atmosphere  (s) 

Range  (in  inertial  system)  reached  during  flight  (m) 

called  for  each  TBM  at  desired  time 


Initial  Threat  East  Longitude  (rad) 

Initial  Threat  North  Latitude  (rad) 

(3  by  1)  Initial  Threat  Velocity  (m/s) 

Initial  Threat  State  Error  Covariance 
Initial  Threat  Altitude  (m) 

Altitude  of  Impact  Area  (m) 

Altitude  where  air  density  is  1/e  of  surface  value  (m) 
Drag  Coefficient  (1/m) 

Time  Increment  for  Numerical  Integration  (s) 
Parameter  Needed  in  Runge  Kutta  Integration 
Parameter  Needed  in  Runge  Kutta  Integration 
Parameter  Needed  in  runge  Kutta  Integration 


Impact  Point  East  Longitude  (rad) 
Impact  Point  North  Latitude  (rad) 
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Ci  (2  by  2)  Covariance  of  Impact  Position  Error  (m?) 

TF  Expected  Flight  Time  (s) 

gentraj  -  called  at  beginning  of  run 

This  is  supplied  by  used  to  generate  a  number  n  of  threat  trajectories.  These  are  used 
to  input  initial  threat  states  and  covariance  matrices  at  different  times. 

Inputs 

none 

Outputs 


n 

Number  of  threat  trajectories 

na 

Number  of  ASCM  Trajectories 

nt 

Number  of  TBM  Trajectories 

nd 

Maximum  number  of  time  steps  in  trajectories 

(nd  by  n  by  6)  Array  of  threat  states 

(nd  by  n  by  6  by  6)  Array  of  covariance  matrices 

Mvec 

(n  by  1)  Threat  Menace  parameter 

kepler  -  a  function  to  find  the  flight  time  and  arc  traversed  in  Keplerian  orbit 

Inputs 

zO 

Initial  Threat  Altitude  (m) 

zi 

Altitude  of  Targeted  Area  (m) 

rOdotcel 

(3  by  1)  Velocity  Vector  in  Celestial  System 

Outputs 

thetakep 

Arc  Distance  Traversed  in  Keplerian  Orbit 

TFkep 

Time  of  Flight  in  Keplerian  Orbit 

makeJ21  -  a  function  that  generates  the  transformation  matrix  J21 

Inputs 

phi  1, LI 

Geographic  Position  of  First  Cartesian  System  (rad) 

phi2,L2 

Geographic  Position  of  Second  Cartesian  System  (rad) 

Output 

J21  (3  by  3)  Transformation  Matrix 
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overlap  -  a  function  that  computes  TBM  impact  probability- 


inputs 

xO,yO  Coordinates  of  Center  of  Defended  Area  (m) 
asemi  Semimajor  Axis  of  Elliptical  Region  (m) 

bsemi  Semiminor  Axis  of  Elliptical  Region  (m) 

phia  direction  of  line  of  apsides  (rad) 

xi,yi  Coordinates  of  Center  of  TBM  Footprint  (m) 

Ci  Covariance  Matrix  of  Footprint  PDF  (m2) 

Output 

prob  Probability  of  Impact  Within  Defended  Area 


simptraj  a  function  called  for  each  (ASCM  -  Defended  Point)  pair  at  desired  time 
Inputs 


rOvec 

rOdotvec 

bm 

phim 

rsvec 

Outputs 

flag 

etarcvec 

etartvec 

etathetat 

eta 

s 


(3  by  1)  Initial  threat  position  (m) 

(3  by  1)  Initial  threat  velocity  (m/s) 
Minimum  allowed  turn  radius  (m) 
Maximum  allowed  heading  error  (rad) 
(3  by  1)  targeted  point  coordinates  (m) 


Indicator  of  ETA  trajectory  validity 
(3  by  1)  Position  of  turn  arc  center  (m) 
(3  by  1)  Position  of  transition  point  (m) 
length  of  turn  arc  (rad) 

ETA  arrival  time  (s) 

Length  of  CV  segment  (m) 


taskfors  -  to  be  written  by  user  and  called  at  desired  times 
Inputs 

none 

Outputs 

m  Number  of  Defended  Assets 

mp  Number  of  Defended  Points 
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ma 

rsvec 

rOvec 

avec 

bvec 

phiavec 

Vavec 

Vpvec 


Number  of  Defended  Areas 

(mp  by  3)  Positions  of  Defended  Points  (m) 

(ma  by  3)  Centers  of  Defended  Areas  (m) 

(ma  by  1)  Semimajor  axes  of  Defended  Areas  (m) 
(ma  by  1)  Semiminor  axes  of  Defended  Areas  (m) 
(ma  by  1)  Directions  of  Lines  of  Apsides  (rad) 

(m  by  1 )  Actual  Values  of  Defended  Assets 
(m  by  1)  Perceived  Values  of  Defended  Assets 


timorder  -  a  function  that  orders  the  ETA  arrival  times 


Inputs 

k 

m 

T(l:n,l:m) 

Va(l:m) 

E(l:n,l:m) 


threat  index 

(or  numit)  Number  of  Defended  Assets 
ETA  Times  of  Arrival 
True  Values  of  Defended  Assets 
Expected  Loss  of  Value 


Outputs 

That(l:n,l:m) 

Vhat(l:n,l:m) 

Ehat(l:n,l:m) 

Etil(l:n,l:m) 


Ordered  ETA  Times  of  Arrival  (s) 
Ordered  True  Values  of  Defended  Assets 
Ordered  Expected  Loss  of  Value 
Maximum  Possible  Value  Loss 


For  computation  of  Maximum  Expected  Value  Loss,  Va*P  is  inputted  instead  of  VA, 
and  the  outputted  Etil  contains  MEVL  values. 


twoarc  -  a  script  file  for  displaying  two-arc  trajectories 


Inputs 

rOvec 

rOdotvec 

rsvec 

nomega 


(3  by  1)  Initial  Threat  Position  (m) 

(3  by  1)  Initial  Threat  Velocity  (m/s) 
(3  by  1)  Targeted  Point  Position  (m) 
(or  na)  Number  of  Turning  Rates 


Outputs 

none 
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uodeint  -  a  function  based  on  the  Runge  Kutta  ODE  solver,  “odeint”  (see  Press  et  al,  1986) 
Inputs 

yinout  (4  by  1)  Vector  of  Dependent  Variables 
tl,t2  Beginning  and  End  Times  of  Stepped  Interval 
epsi  Parameter  Governing  Integration  Process 

hl,hmin  Parameters  Governing  Integration  Process 

Outputs 

yinout  (4  by  1)  Vector  of  Dependent  Variables 
nok,nbad  Parameters  Displaying  Result  of  Step 

urkqc  -  a  function  used  in  Runge  Kutta  process,  based  on  “rkqc”  (see  Press  et  al,  1986) 

Inputs 

yin  (4  by  1)  Vector  of  Dependent  Variables 

tin  (4  by  1)  Independent  Variable 

delta  Modified  Drag  Coefficient 

hatm  Altitude  where  Air  Density  is  1/e  of  Surface  Value 
epsi,yscal  Parameters  Governing  Stepping  Process 

Outputs 

yout  Vector  of  Dependent  Var  iables 

tout  Independent  Variable 

hdid  Parameter  Displaying  Result  of  Step 

hnext  Parameter  Displaying  Result  of  Step 

urk4  -  function  used  in  Runge  Kutta  process,  based  on  “rk4”  (see  Press  et  al,  1986) 

Inputs 

yin  (4  by  1)  Vector  of  Dependent  Variables 

tin, tout  Times  at  Beginning  and  End  of  Step 

delta  Modified  Drag  Coefficient 

hatm  Altitude  where  Air  Density  is  1/e  of  surface  Value 

Outputs 

yout  Vector  of  Dependent  Variables 
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Figure  10.1.  Organization  of  the  Computations  for  Threat  Prioritization. 
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CHAPTER  11 

SUMMARY  AND  CONCLUSIONS 

This  report  considers  the  case  where  a  number  of  airborne  threats  may  be  targeting  a 
number  of  defended  assets.  The  threats  are  assumed  to  be  one  of  two  types:  Anti-Ship 
Cruise  Missiles  (ASCMs)  or  Theater  Ballistic  Missiles  (TBMs).  The  defended  assets  are  of 
two  types:  points  (e.g.,  ships)  or  areas  (e.g.,  cities).  Each  threat  is  evaluated  to  provide 
parameters  to  be  passed  to  the  Weapons  Control  System  to  aid  in  prioritization  and  en¬ 
gagement.  The  threat  evaluation  is  usually  based  on  the  latest  state  estimate  and  whatever 
classification  information  is  available. 

ASCMs  are  treated  by  constructing  a  simplified  constant-speed  trajectory  to  each  de¬ 
fended  asset  to  estimate  the  Earliest  Time  of  Arrival  (ETA).  These  trajectories  obey  con¬ 
straints  limiting  the  lateral  acceleration  magnitude  and  heading  error.  Objective  Probabil¬ 
ities  are  estimated  for  each  ASCM-asset  pairing,  based  on  the  asset’s  value,  ETA,  and  the 
energy  requirement  for  reaching  the  asset. 

TBMs  are  assumed  unable  to  maneuver.  Uncertainty  in  the  impact  point  arises  because 
of  state  estimation  errors  and  lack  of  knowledge  of  the  physical  characteristics  of  the  vehicle. 
This  uncertainty  is  estimated  to  produce  an  impact  Probability  Density  Function  (PDF), 
or  footprint  ,  representing  the  place  where  the  TBM  might  land.  Objective  Probabilities 
for  each  TBM-asset  pairing  are  estimated  on  the  basis  of  the  asset’s  value  and  the  footprint 
PDF  at  the  asset’s  position. 

The  evaluation  of  each  threat  is  expressed  in  terms  of  two  functions  of  time:  the  Max¬ 
imum  Expected  Value  Loss  (MEVL)  and  the  Maximum  Possible  Value  Loss  (MPVL),  The 
Objective  Probabilities,  asset  values,  ETAs,  and  other  information  are  used  to  calculate  the 
MEVL  and  the  MPVL.  At  any  time  the  MEVL  and  the  MPVL  for  a  given  threat  express  in 
slightly  different  ways  the  value  loss  that  would  be  incurred  if  the  threat  were  to  leak. 

It  is  recognized  that  the  assumptions  upon  which  the  hypothesized  predictive  trajectories 
are  based  may  not  be  universally  applicable.  The  main  assumptions  have  been  (a)  constant 
speed,  (b)  a  limited  HE,  and  (c)  a  limited  “g”  capability.  Constant  speed  results  from 
neglect  of  transfer  between  potential  and  kinetic  energy  upon  descent  or  ascent,  as  well  as 
from  disregard  of  drag  losses  and  the  possibility  of  powered  flight.  A  limited  HE  assumes  that 
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the  threat  must  continually  be  in  a  homing  (seeking)  mode.  It  is  likely  that  less  restrictive 
predictive  models  will  provide  better  estimates  of  ETA  and  a  more  useful  threat  prioritization 
algorithm. 
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APPENDIX  A 

ESTIMATION  OF  TBM  DRAG 
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A.l  Introduction 

Air  has  a  negligible  effect  on  a  Theater  Ballistic  Missile  (TBM)  during  most  of  its  flight. 
Only  the  initial  and  final  phases  of  flight  are  influenced  by  aerodynamic  effects.  During  the 
initial  phase  the  TBM  is  usually  subject  to  a  large  thrust  force  that  makes  it  difficult  to 
estimate  the  aerodynamic  effects.  Therefore,  the  concepts  discussed  here  apply  only  during 
the  final  seconds  of  flight  when  it  may  be  too  late  to  take  defensive  action. 

For  aerodynamic  forces  on  a  body,  the  force  is  a  function  of  the  speed.  If  the  force  is  a 
function  of  the  speed  alone,  the  force  law  can  be  expressed  as 


n  rn/  %  rfW 

Fq  +  F(v )  =  m— 
at 


(Al) 


Fq  +  F(v)  =  mv- 


(A- 2) 


Fo  is  any  constant  force  that  doesn’t  depend  on  v.  For  fluid  resistance  a  fair  approximation 
to  F(v)  is 

F(v )  =  -c\v  —  C2v\v\  (A3) 

where  ci  and  C2  are  constants  that  depend  on  the  size  and  shape  of  the  object.  For  spherical 
objects,  the  approximate  values  of  the  constants  in  (A.3)  are  c\  =  c\D  and  C2  =  b^D2  where 
D  is  the  diameter  of  the  sphere  and  Ci  =  1.55  X  10~4  kg/sec  and  £2  =  .22  kg/m.  For  two 
cases,  there  is  a  limiting  form  of  the  drag  problem  depending  on  whether  the  speed  is  large 
or  small.  For  small  v  the  limiting  form  is 


F(v)  =  kpv  —  mgz 

(z  is  the  upward  unit  vector)  while  for  large  v,  the  limiting  form  is 

F(v)  =  kp  |F|  v  —  mgz 


The  third  alternative  for  the  friction  law  is 


(A- 4) 


(A5) 


F(v)  =  kpv 2 

for  the  x  and  y  components,  while  the  z  component  is 

F(v)  —  kpv 2  —  mg 


(A.6) 


(A.  7) 
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Closed  form  solutions  can  be  obtained  for  (A.4)  and  for  (A.6-A.7)  but  not  for  (A. 5).  An 
alternative  to  (A.4)  that  is  sometimes  used  is  to  multiply  v  by  e~z!h  in  (A.4).  What  we  are 
interested  in  determining  is  estimators  of  the  ballistic  coefficient.  For  the  dynamics  laws  that 
have  closed  form  solutions  to  the  equations  of  motion  can  be  used  as  a  means  for  estimating 
the  ballistic  coefficient.  For  other  dynamics,  alternative  estimators  are  derived. 


A. 2  Solution  of  Equations  of  Motion 

For  motion  along  the  2  axis,  in  the  quadratic  case  there  are  two  differential  equations  to 
consider:  upward  motion  (positive  v,  minus  sign)  and  downward  motion  (negative  v,  plus 
sign) 

1  2  < 

—mg  ±  c2v  =  m- 


dt 


The  solution  to  (A. 8)  is 

fv  mdv  (  (v0\  fvW  , 

t  —  to  =  I  - r  =  r  (  arctan  I  —  I  —  arctan  I  —  I  1  ( ris ) 

Jv0-m9-C2V*  \  \vtJ  \vt/J 


-*-f  z; 

Jvi]  ' 


m  dv 


=  r 


arctanh 


V 


(£)-arctanh(£))  (fal) 


lv 0  - mg  +  c2v 2 

where  r  =  ■\fmjc2g  and  v*  =  \Jrng]~c2.  Solving  these  equations  for  v  gives 


v  =  v<tan 


v  =  — ujtanh 


to  -t  vo 

- h  arctan — 

r  Vi 

t  -  t'r 


-  —  arctanh— 
r  vt 


(ris) 

(fal) 


The  position  as  a  function  of  time  is  therefore 

'to  -  t 


Z  —  ZO  =  VfT log  (  COS 


UO 


+  arctan —  )  1  (?-is) 

T  Vf 


and 


(A.S) 

(A.9) 

(4.10) 


(4.11) 
C 4.12 ) 

(-4-13) 


z  —  zq  —  v^rlog  ^cosh  — -  —  arctanh— ^  (/a/)  (^4. 14) 

Note  for  both  of  these  cases,  the  x  and  y  positions  are  the  same  as  those  in  (A.11-A.12). 
Expanding  (A.13-A.14)  using  the  standard  identities  gives 


(z  -  zo)/(vtT)  =  -log^^  +  log  vtcos  ~  vosin 

(z  -  zo)/(vtT)  =  -log ^2  -  vl  +  log  Dfcosh  _  vosinh  ^ 


(A.15) 
(A.  16) 
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For  the  x  and  y  positions,  the  equations  of  motion  follow  from  solving  (A.5),  which  gives 


x  —  xq  =  —  In  (bxvjism(a)t  +  1) 


(AA7) 


V-yo  =  r  In  (byVRCOs(a)t  +  1) 


(A18) 


We  will  now  assume  both  to  and  are  both  zero  for  the  sake  of  simplicity.  The  second  order 
expansions  for  rising  and  falling  objects,  (A.15-A.16),  become 


z  -  ZQ  =  u^rlog 


,vt  +  vl 


+  v„t+LAjL_'Z 

T  \Vt  2 


(A.  19) 


Z  —  ZQ  =  VfTlog 


'  Vi  ~  Vl 


+  u0t  +  vt  +  2—  ( fal ) 

IT  Vt 


(A.20) 


Thus,  the  effect  of  a  quadratic  drag  term  is  to  induce  an  ’acceleration’  whose  magnitude 
using  the  relationship  vt/r  =  g,  becomes  either  2vq/ gr2  ±  g.  Approximating  the  solutions 
for  the  x  and  y  components  gives 


x  -  z°  =  [uR&ssin^jf  -  \v2Rb2xsm2{a)t2  +  ...] 


(A.21) 


y-yo  =  Y  [u/j6j,cos(q!)  -  \v2Rb2 cos2 (o:)f2  +  ...] 
Note  if  we  retain  first  order  terms  only,  the  slant  range  is 

RH{t)  =  J{x-  *o)2  +  (y  -  m)2  =  vRt 


(A.  22) 


(A.23) 


Thus,  under  the  circumstance  that  the  drag  is  small,  one  could  ignore  the  effects  of  drag  on 
the  slant  range  velocity.  Thus  it  suffices  to  design  a  two  state  filter  to  handle  the  slant  range 
when  this  approximation  is  true,  A  more  realistic  assumption  is  to  keep  the  second  order 
term  because  b  is  in  the  range  50-500.  For  this  situation,  with  the  assumption  bx  =  by ,  it  is 
better  to  go  to  third  order  filter. 

A. 3  Ballistic  Estimators 

These  formulas  enable  one  to  determine  the  ballistic  coefficient  from  tracking  data,  they 
do  not  work  for  the  types  of  drag  functions  that  do  not  have  analytic  solutions.  The  method 
is  general  enough  to  work  for  any  dispersive  system  of  the  type  we  are  discussing.  An 
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estimator  can  be  derived  based  on  the  following  observation:  Any  dispersive  system  that 
can  be  written  as 

7t{T  +  v)=7tEc  =  F'v  (A24) 

where  F’  is  the  portion  not  derivable  from  a  potential,  T  -f  V  is  kinetic  plus  potential  energy 
or  the  total  conservative  energy  Ec,  and  F!  is  the  portion  of  the  force  not  derivable  from  a 
potential  or  the  power  delivered  by  an  effective  force  F' . 

For  our  application,  the  normalized  Ec  is 

Ec  =  \ v 2  +  gz  (A.25) 

There  are  several  possibilities  for  the  effective  force,,  depending  on  the  circumstances.  An  es¬ 
timator  can  be  determined  independent  of  the  details  of  what  F'(v)  functional  form  provided 
it  is  assumed  that 

F'(v)  =  kpf(v)  (A.  26) 

Substituting  (A.25)  and  (A. 26)  into  (A. 24)  gives 

^(±v2  +  gz)  =  kpvf(v) 
or 


d(\v2  A  gz)  =  kpvf{v)dt  (A. 27) 

To  derive  an  estimator  from  (A. 27),  the  derivatives  and  the  variables  are  made  discrete. 

d  —*  A{  =  (Si  —  5j-i) 


so  (A. 27)  becomes 

Ei  -  Ei-i  =  kipiVif(vi)(U  -  ti- 1) 

(A. 28) 

Summing  over  i  gives 

n 

En  -  E0  =  ^  bpiVif(vi)(ti  -  ti- 1) 

(A.29) 

t=i 


If  one  defines  D  as  the  negative  of  the  right  side  of  (A. 29),  and  recognizes  that  the  variables 
are  measurements,  the  energy  loss  estimator  LE  is  formed 

LE  =  D  +  En  -  E0  (A. 30) 

which  can  be  used  as  ones’  test  statistic  for  a  hypothesis  test.  General  Electric  (Reifler)  has 
attempted  this.  There  are  other  alternatives  to  the  energy  estimator  in  (A.35)  besides  those 
already  mentioned. 
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Another  method  is  to  solve  the  equation  of  motion  for  the  drag  coefficient,  obtaining 


8  =  v  2ezlh 


dv 

dt 


g  sm  7 


(A.31) 


The  change  in  threat  speed  An  is  observed  a  number  n  of  times  along  the  TBM  trajectory, 
resulting  in  a  set  of  values 

Avj,  Atj,  J  =  1, n  (A.  32) 

At  each  observation  an  estimate  of  the  drag  coefficient 


8  +  v-  2ezi/h 


Avgj 

~At~ 


(A.  33) 


is  obtained,  where  8  is  the  true  value  and  A vej  are  the  errors  in  the  observed  speed  increments 
with  standard  error  cr&v.  The  weighted  mean  of  these  estimates  is 


8  = 


J2  Sjvje  zilhAtj 


I< 


=  8- 


£At>, 


ej 


I< 


k  =  'E 


(A.34) 


This  weighting  provides  an  unbiased  estimate  of  the  drag  coefficient.  To  find  the  standard 
error  of  this  estimate  it  is  necessary  to  find  the  standard  error  in  the  speed  estimate.  This  is 
cumbersome  in  general  even  when  the  error  in  the  estimates  of  the  threat  state  are  assumed 
Gaussian.  However,  usually  the  threat  speed  greatly  exceeds  its  error.  Accordingly,  it  is 
assumed  that 

|t»-  (i>}|  <  H>l  (A.35) 

where  the  angle  brackets  denote  expected  value.  The  threat  speed  v  =  |u|  is  expressed  by 


v  = 


(A.36) 


where  the  magnitude  of  the  vector  within  the  parentheses  is  assumed  small.  To  first  order 
in  the  small  quantity  this  becomes 


»  =  W0+2W -2W(<?-(ff))  =  (A37) 

and 

(v)  =  l(u)I  (A. 38) 

Thus,  since  v  is  a  linear  combination  of  the  components  of  v,  v  is  Gaussian  if  v  is  Gaussian. 
The  error  in  the  speed  estimates  is 

av  =  <[«  “  (v)]2)  =  (A.39) 
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where  C  is  the  error  covariance  of  the  velocity  estimates, 

c=(|tr-(S)][0-<v»T)  (-4.40) 

Estimates  of  the  speed  increments  are  not  independent  even  when  the  speed  errors  are 
independent.  If  one  increment  is  overestimated,  the  following  increment  is  likely  to  be 
underestimated.  However,  the  error  in  a  speed  increment  is  independent  of  the  errors  in  all 
the  non-contiguous  speed  increments.  The  correlation  between  contiguous  increment  errors 
has  little  effect,  and  the  relation  between  them  and  those  of  the  speed  estimates  is 


=  2  <rj 


(-4.41) 


The  error  in  the  drag  estimate  is  then 


°s  =  ([<*-  <<*>]2)  = 


EAuej-12\  _2 QfCJft 

K  J  / 


(A.  42) 


The  foregoing  analysis  is  based  on  a  number  of  assumptions,  some  of  which  have  been 
discussed  previously.  One  is  that  the  air  density  is  a  simple  exponential  function  of  altitude. 
A  more  refined  approach  would  be  to  take  the  actual  density  profile  existent  at  the  place  and 
time.  Another  is  that  the  speed  errors  are  independent  of  time  (and  altitude),  whereas  in 
actuality  these  errors  are  likely  to  increase  as  the  air  density  and  drag  increase.  Still  another 
is  that  the  expected  velocity  varies  little  during  the  time  the  observations  a,re  made. 
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This  scenario  illustrates  the  time-ordered  threat  evaluation  parameters,  the  Maximum 
Possible  Value  Loss  (MPVL)  and  the  Maximum  Expected  Value  Loss  (MEVL),  defined  in 
Section  9.5  of  the  main  text.  There  is  a  task  force  consisting  of  fifteen  defended  points 
and  four  Anti-Ship  Cruise  Missile  (ASCM)  trajectories,  as  shown  in  Figure  B.l.  The  size 
of  each  dot  representing  a  defended  point  is  proportional  to  its  value.  Here,  there  is  no 
distinction  between  actual  and  perceived  value.  The  ASCM  trajectories  are  numbered  at 
the  point  representing  their  initial  position.  The  trajectories  are  shown  in  plan  view,  the 
threat  positions  being  projected  vertically  onto  the  horizontal  plane. 

The  evolution  of  the  Earliest  Times  of  Arrival  (ETAs)  is  shown  in  Figure  B.2  for  each 
of  the  four  threats.  There  are  fifteen  curves,  one  for  each  of  the  defended  points.  The 
curve  for  each  defended  point  terminates  after  the  threat  is  no  longer  able  to  reach  it.  These 
curves  generally  decrease  with  increasing  actual  time,  but  there  is  sometimes  a  rather  abrupt 
increase  in  arrival  time  when  the  threat  progresses  into  a  position  less  favorable  for  reaching 
that  target. 

The  evolution  of  ETA  drag  loss  is  shown  in  Figure  B.3  for  each  of  the  four  threats.  These 
curves  generally  display  the  same  nature  as  those  for  the  arrival  times.  The  rate  of  drag  loss 
increases  with  decreasing  radius  of  curvature,  and  thus  the  curves  rise  whenever  the  threat 
needs  to  pull  a  higher  lateral  acceleration  to  reach  the  targeted  point. 

The  evolution  of  the  Objective  Probabilities  for  Threat  =^1  is  shown  in  Figures  B.4a-e 
at  five  different  times.  The  size  of  the  dots  representing  the  defended  points  is  proportional 
to  the  Objective  Probability.  Points  no  longer  reachable  by  the  threat  are  indicated  by 
their  absence  from  the  figure.  The  Objective  Probabilities  for  Threat  #2  are  shown  at  three 
different  times  in  Figure  B.5a-c,  and  for  Threats  #3  and  #4  in  Figures  B,6a-b  and  B.7a-c. 

The  MPVL  is  displayed  for  the  four  threats  at  three  different  times  in  Figures  B.8a-c. 
The  MEVL  is  similarly  displayed  in  Figures  B.9a-c.  The  amount  of  information  represented 
in  these  displays  is  small,  as  only  the  coordinates  of  the  corner  points  of  the  staircase  repre- 
sentation  are  needed. 

Figure  B.10  displays  the  comprehensive  MPVL  and  MEVL  functions  for  all  the  threats. 
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CO 


B-4 


Figure  B.l.  An  example  scenario  of  fifteen  defended  points  and  four  ASCM  trajectories  in 
plan  view.  The  trajectories  are  projected  vertically  onto  the  sea  surface.  The  size  of  the 
dots  representing  the  defended  points  is  proportional  to  their  actual  values. 
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Figure  B.4a.  Objective  Probabilities  for  the  fifteen  defended  points,  indicated  by  size  of  the 
dots,  for  threat  $1  at  time  7  sec. 


for  time  =  8,  traj  1 
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Figure  B.4b.  Objective  Probabilities  for  the  fifteen  defended  points,  indicated  by  size  of  the 
dots,  for  threat  #1  at  time  8  sec. 


for  lime  =  9,  Iraj  1 
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Figure  B.4c.  Objective  Probabilities  for  the  fifteen  defended  points,  indicated  by  size  of  the 
dots,  for  threat  #1  at  time  9  sec. 


for  lime  =  10,  traj  1 
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Figure  B.4d.  Objective  Probabilities  for  the  fifteen  defended  points,  indicated  by  size  of  the 
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for  lime  *  7,  traj  2 
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Figure  B.5a.  Objective  Probabilities  for  the  fifteen  defended  points,  indicated  by  size  of  the 
dots,  for  threat  ,#2  at  time  7  sec. 
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Figure  B.5b.  Objective  Probabilities  for  the  fifteen  defended  points,  indicated  by  size  of  the 
dots,  for  threat  ,#2  at  time  12  sec. 


for  lime  =  13,  traj  2 
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Figure  B.5c.  Objective  Probabilities  for  the  fifteen  defended  points,  indicated  by  size  of  the 
dots,  for  threat,  #2  at  time  13  sec. 
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Figure  B.6a.  Objective  Probabilities  for  the  fifteen  defended  points,  indicated  by  size  of  the 
dots,  for  threat  #3  at  time  13  sec. 


for  lime  *19,  traj  3 


NSWCDD/TR-00/46 


B-16 


Figure  B.6b.  Objective  Probabilities  for  the  fifteen  defended  points,  indicated  by  size  of  the 
dots,  for  threat,#3  at  time  19  sec. 


for  time  =  4,  traj  4 
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Figure  B.7a.  Objective  Probabilities  for  the  fifteen  defended  points,  indicated  by  size  of  the 
dots,  for  threat,  #4  at  time  4  sec. 


for  time  =  14,  traj  4 
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Figure  B.7b.  Objective  Probabilities  for  the  fifteen  defended  points,  indicated  by  size  of  the 
dots,  for  threat,  #4  at  time  14  sec. 


for  times  16Jraj4 
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Figure  B.7c.  Objective  Probabilities  for  the  fifteen  defended  points,  indicated  by  size  of  the 
dots,  for  threat.  #4  at  time  16  sec. 
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Figure  B.8a.  The  Maximum  Possible  Value  Loss  (MPVL)  for  the  four  threats  of  the  scenario 
of  Figure  B.l  at  time  t  =  0. 
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Figure  B.8b.  The  Maximum  Possible  Value  Loss  (MPVL)  for  the  four  threats  of  the  scenario 
of  Figure  B.l  at  time  t  =  1. 
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Figure  B.8c.  The  Maximum  Possible  Value  Loss  (MPVL)  for  the  four  threats  of  the  scenario 
of  Figure  B.l  at  time  t  =  2. 
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Figure  B.9a.  The  Maximum  Expected  Value  Loss  (MEVL)  for  the  four  threats  of  the 
of  Figure  B.l  at  time  t=0. 
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Figure  B.9b.  The  Maximum  Expected  Value  Loss  (MEVL)  for  the  four  threats  of  the  scenario 
of  Figure  B.l  at  time  t=l. 


MAXIMUM  EXPECTED  VALUE 
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Figure  B.9c.  The  Maximum  Expected  Value  Loss  (MEVL)  for  the  four  threats  of  the  scenario 
of  Figure  B.l  at  time  t=2. 
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Figure  B.10.  Evolution  in  time  of  the  Comprehensive  MPVL  and  MEVL. 
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APPENDIX  C 

VARIABLE  NAMES 
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I 


Comprehensive  List  of  Variable  Names 


The  units  indicated  apply  only  to  the  variable  name  in  the  computer  program. 

Chapter  2 


3 

Va 

actual  value  of  a  defended  asset  (ndim) 

l/.W 

3 

vP 

perceived  value  of  a  defended  asset  (ndim) 

$ 

elong 

East  longitude  (rad) 

L 

nlat 

North  latitude  (rad) 

r 

r 

radial  coordinate  (m) 

<fa 

elongO 

East  longitude  of  point  of  tangency  (rad) 

Lo 

nlatO 

North  latitude  of  point  of  tangency  (rad) 

u 

u 

East  velocity  component  (m/s) 

V 

V 

North  velocity  component  (m/s) 

w 

w 

upward  velocity  component  (m/s) 

XI 

Cartesian  component  in  East  direction  (m) 

yi 

Cartesian  component  in  North  direction  (m) 

Z\ 

Cartesian  component  in  upward  direction  (m) 

<pl 

East  longitude  at  tangent  point  of  Local  Coordinate 
System  (LCS)  (rad) 

Li 

L  at  tangent  point  of  LCS  (rad) 

a 

a 

Earth  radius  (m) 

U\ 

ul 

velocity  component  in  xl  direction  (m/s) 

Vl 

vl 

velocity  component  in  yl  direction  (m/s) 

W\ 

wl 

velocity  component  in  zl  direction  (m/s) 

X\ 

xtil 

Cartesian  coordinate  in  second  system  (m) 

m 

ytil 

Cartesian  coordinate  in  second  system  (m) 

h 

ztil 

Cartesian  coordinate  in  second  system  (m) 

Xi 

xhat 

Cartesian  coordinate  in  third  system  (m) 

m 

yhat 

Cartesian  coordinate  in  third  system  (m) 

z\ 

zhat 

Cartesian  coordinate  in  third  system  (m) 
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(Note  that  (x2,y2,Z2)  and  (x2,y2,Z2)  are  similarly  defined  in  another  LCS  at  another 
point  of  tangency.) 


Jjk 

Jacobian  matrix  (var) 

“9” 

e 

suffix  indicating  geographic  coordinates 

a 

asemi 

semimajor  ellipse  axis  (m) 

b 

bsemi 

semiminor  ellipse  axis  (m) 

e 

e 

eccentricity  (ndim) 

(zo,yo) 

(xO,yO) 

LCS  location  of  ellipse  center  (m) 

(xfhVfi) 

(xfl,yfl) 

LCS  location  of  ellipse  focus  (m) 

(z/2,2//2) 

(xfl,yfl) 

LCS  location  of  other  focus  (m) 

0a 

phia 

trend  (direction)  of  line  of  apsides  (rad) 

e 

theta 

angular  polar  coordinate  (rad) 

j 

j 

index  for  defended  asset 

k 

k 

threat  index 

3 

Wa 

actual  value  per  unit  area  (m-2) 

< 

Wp 

perceived  value  per  unit  area  (m-2) 

Chapters  3  and  4 

n 

n 

number  of  threats  (ndim) 

na 

na 

number  of  cruise-missile  tracks  (ndim) 

nt 

nt 

number  of  ballistic  missile  tracks  (ndim) 

m 

m 

number  of  defended  assets  (ndim) 

rrip 

mp 

number  of  defended  points  (ndim) 

ma 

ma 

number  of  defended  areas  (ndim) 

nd 

nd 

number  of  data  points  in  trajectory  (ndim) 

fm 

fm 

maneuver  factor  (ndim) 

ro 

rOvec 

initial  threat  position  (m) 

— # 

ro 

rOdotvec 

initial  threat  velocity  (m/s) 

rs 

rsvec 

asset  position  (m) 

Vt 

Vt 

threat  speed  (m/s) 

bm 

bm 

minimum  allowed  turn  radius  (m) 

% 

am 

maximum  allowed  lateral  acceleration  (m/s2) 

0 

phi 

heading  error  (rad) 

0m 

phim 

maximum  allowed  heading  error  (rad) 

Mk 

M(k) 

threat  “Menace”  (ndim) 
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A,B 

A,B 

undetermined  constants 

rc 

rcvec 

center  of  curvature  (m) 

E,F 

E,F 

undetermined  constants 

n 

rtvec 

transition  point  (m) 

A 

Ahat 

undetermined  constant 

B 

Bhat 

undetermined  constant 

C 

Chat 

undetermined  constant 

Bt 

thetat 

arc  length  traversed  (rad) 

6 

theta 

arc  length  (rad) 

b 

b 

radius  of  curvature  (m) 

Os 

thetas 

another  arc  length  (rad) 

n 

rivec 

splash  point  or  impact  point  (m) 

a,  b 

a,b 

legs  of  a  right  triangle  (m) 

s 

s 

length  of  a  CV  segment  (m) 

Chapter  6 

fd  fd 
fl  fl 

k  kappa 
Cm  Cdo 
K  K 

p  rho 

S  S 

m  m 
ro  rddotvec 

Fd 

Fl 

Th 

9 

k 

u 

a  a 
a 

L0  Lo 


drag  coefficient  (ndim) 

lift  coefficient  (ndim) 

a  physical  parameter  of  air  (km/m) 

drag  coefficient  at  no  lift  (ndim) 

an  aerodynamic  parameter  for  a  threat  (ndim) 

air  density  (kg/m3) 

threat  cross  section  area  (m2) 

threat  mass  (kg) 

initial  threat  acceleration  (m/s2) 

drag  force  (kg  m/s2) 

lift  force  (kg  m/s2) 

thrust  force  (kg  m/s2) 

acceleration  of  gravity  (m/s2) 

upward  unit  vector  (ndim) 

unit  vector  in  direction  of  lift  (ndim) 

acceleration  magnitude  (m/s2) 

acceleration  (m/s2) 

drag  loss  (kg  m2/s2) 
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Chapter  7 

t 

t 

time  (s) 

n 

omega 

Earth  rotation  rate  (1/s) 

Oc 

thetac 

arc  length  traversed  (rad) 

T 

T 

flight  time  (s) 

Xo 

XO 

flight  range  in  Keplerian  trajectory  (m) 

Xa 

Xa 

flight  range  in  air  (m) 

To 

TO 

flight  time  in  Keplerian  trajectory  (s) 

Ta 

Ta 

flight  time  in  air  (s) 

Vo 

VO 

initial  TBM  speed  (m/s) 

Re 

radial  coordinate  of  defended  area  (m) 

X 

lambda 

parameter  in  Keplerian  calculation  (ndim) 

V 

mu 

parameter  in  Keplerian  calculation  (m3s-2) 

h 

hatm 

altitude  where  air  density  is  1/e  surface  value  (m) 

cr 

sigma 

parameter  in  Keplerian  calculation  (ndim) 

Tf 

TF 

Keplerian  flight  time  (s) 

7 

gamq 

flight  path  angle  (rad) 

70 

gamqO 

initial  flight  path  angle  (rad) 

<5 

Gray’s  drag  parameter  (1/m) 

C 

zeta 

square  of  threat  speed  (m2/s2) 

a 

alpha 

an  angle  in  spherical  triangle  (rad) 

Xfall 

horizontal  range  of  fall  through  air  (m) 

Tfall 

flight  time  of  fall  through  air  (s) 

espi 

parameter  in  Runge  Kutta  integration 

hi 

parameter  in  Runge  Kutta  integration 

hmin 

parameter  in  Runge  Kutta  integration 

Chapter  8 

J 

J 

Jacobian  (matrix) 

CO 

CO 

threat-state  error  covariance  matrix 

Ci 

Ci 

TBM  footprint  covariance  matrix  (m2) 

6 

deltab 

expected  drag  parameter  value  (1/m) 

<?6 

sigmad 

standard  error  in  drag  parameter  (1/m) 

Tfi 

TFI 

flight  time  with  different  drag  (s) 

Q\ 

arc  length  with  different  drag  (rad) 

C-6 

Cs  Cs 
Cd  Cd 

Chapter  9 

Pkj  P(kj) 
Tki  T(kj) 
Lkj  L(k,j) 
Ekj  E(kj) 
Eg’  Ep(kj) 
Pkj 

i  iunit 


Conflicts 


a 

b 
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footprint  covariance  from  state  errors  (m^) 
footprint  covariance  from  drag  error  (m2) 


Objective  Probability  (ndim) 

ETA  arrival  time  (s) 

ETA  drag  loss  (kg  m2  s-2) 

Maximum  Expected  Value  Loss  (ndim) 
Maximum  Possible  Value  Loss  (ndim) 

Objective  Probability  for  defended  point  (ndim) 
ordered  index  for  defended  asset 


Earth  radius,  semimajor  axis,  acceleration,  leg  of  triangle 
radius  of  curvature,  semiminor  axis,  leg  of  triangle 
heading  error,  E  longitude 
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